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I.  INTRODUCTION 


The  effects  of  thermal  radiation  in  compressible  flow  have 
received  considerable  attention  in  recent  years  as  a  result  of  the 
increasing  flight  speeds  connected  with  space  flight.  Little  work  has 
been  done,  however,  in  multi-dimensional  flow  problems  with  noneauilb- 
rium  heat  transfer  due  to  radiation.  In  fact,  no  solutions  valid  for 
the  complete  range  of  absorption  coefficient  have  been  obtained  for 
such  problems.  The  difficulty  is  due  essentially  to  the  fact  that  the 
governing  equations  are  of  a  complicated  in tegro- differential  form. 

The  present  work  represents  an  attempt  to  understand  various 
aspects  of  multi- dimensional  radiating  gas  flow  with  no  restriction  as 
to  the  magnitude  of  the  absorption  coefficient.  Chapters  II  end  III 
provide  background  material  that  is  needed  for  subsequent  discussion. 

In  particular.  Chapter  II  reviews  briefly  the  fundamentals  of  radiation 
transport.  In  this,  emphasis  is  placed  on  the  particle  characteristics 
of  radiation.  This  approach  enables  us  to  discuss  the  problems  of 
radiation  transport  in  a  manner  analogous  to  that  used  in  kinetic  theory 
of  gases  and  in  neutron- transport  theory. 

The  basic  integro- differential  equations  of  radiation  gas- 
dynamics  are  reviewed  in  Chapter  III.  For  one- dimensional  problems,  the 
integral  term  in  the  differential  equations  can  be  removed  by  the  so- 
called  "exponential  approximation"  suggested  by  Vincenti  and  Baldwin  (1962). 
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As  a  result,  analytical  solutions  valid  for  the  full  range  of  absorption 
coefficient  are  possible  in  the  one-dimensional  case.  Since  the  expo- 

'  4 

nential  approximation  is  related  to  the  spherical-harmonic  approximation 
that  we  shall  discuss  at  length  in  subsequent  chapters,  a  brief  treatment 
of  this  approximation  is  included. 

Although  the  exponential  approximation  does  lead  to  analytical 
solutions  valid  for  the  complete  range  of  absorption  coefficient  in  the 
case  of  one-dimensional  flow,  this  technique  is  not  applicable  to  multi¬ 
dimensional  problems.  We  therefore  look  for  other  approximate  methods 
that  can  be  used  in  multi- dimens iona1  flow.  One  such  method,  known  as 
the  moment  method,  has  been  discussed  recently  by  Traugott  (1963)  and 
by  the  present  author  (1964).  With  this  approximation,  the  exact  radia¬ 
tion  transport  equation  can  be  replaced  by  a  set  of  its  moment  equations. 
This  set  of  approximate  transport  equations,  together  with  the  gas- 
dynamic  equations,  constitute  a  determinate  set  of  purely  differential 
equations  valid  for  the  complete  range  of  absorption  coefficient.  In 
Section  4.1,  the  same  approximate  transport  equations  for  a  multi¬ 
dimensional  radiating  gas  are  rederived  as  the  first  approximation  in 
a  systematic  spherical-harmonic  method  similar  to  that  used  in  astro¬ 
physics  and  neutron- transport  theory. 

In  the  approximate  treatment  a  question  arises  as  to  how  the 
boundary  conditions  on  the  radiation  field  are  to  be  approximated  in  a 
manner  consistent  with  the  spherical- harmonic  approach.  To  handle  this 
question  the  so-called  "Mark"  and ’Barshak"  boundary  conditions  of  neutron- 
transport  theory  are  applied  to  the  present  problem.  As  a  result, 
expressions  for  the  approximate  boundary  condition  consistent  with  the 
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spherical-harmonic  method  are  given  for  a  black  surface  in  Section  4.2. 
With  the  spherical-harmonic  approximation,  the  governing  equations  and 
the  boundary  conditions  are  in  purely  differential  form.  We  therefore 
refer  to  the  present  formulation  as  the  "differential  approximation" 
to  radiating  gas  flow. 

The  relation  of  the  exponential  approximation  and  the  present 
differential  approximation  is  discussed  in  Section  4.3.  It  is  shovn 
that  the  exponential  approximation  with  suitable  choice  of  constants 
in  one  exponential  function  is  completely  equivalent  to  the  differential 
approximation  when  specialized  to  one- dimensional  situations  and  employed 
with  the  Mark  boundary  condition. 

Even  though  the  governing  equations  in  the  present  formulation 
are  in  purely  differential  form,  they  are  nonlinear,  and  analytical 
solution  is  still  difficult.  We  therefore  examine  the  corresponding 
linearized  equations  in  Chapter  V.  As  a  result,  a  single,  linear,  fifth- 
order  partial  differential  equation  in  terms  of  the  velocity  potential, 
with  coordinate  system  fixed  in  the  undisturbed  fluid,  is  obtained.  The 
corresponding  equation  for  a  coordinate  system  fixed  in  a  body  moving 
with  constant  speed  is  found  by  means  of  a  Galilean  transformation. 

These  two  equations  play  the  same  role  as  the  acoustic  equation  and  the 
Prandtl-Glauert  equation  respectively  in  classical  gasdynamics.  Within 
the  framework  of  the  linearized  theory,  an  expression  for  the  boundary 
condition  on  the  wall  temperature  is  given.  This  expression  shows  that, 
except  in  the  limiting  case  of  an  opaque  gas,  a  temperature  Jump  exists 
at  the  wall. 
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To  illustrate  the  present  formulation  of  radiating  gas  flow, 
several  simple  examples  are  considered  in  Chapter  VI.  In  particular, 
the  propagation  of  linearized  one- dimensional  waves,  solved  previously 
with  the  exponential  approximation,  is  reconsidered  on  the  basis  of  the 
differential  approximation.  It  is  shown  that  the  results  obtained  from 
the  two  methods  of  approximation  are  essentially  equivalent.  Finally, 
with  the  desire  to  provide  a  simple  analytical  solution  for  multi¬ 
dimensional  radiating  flow,  the  two-dimensional  problem  of  flow  over 
a  wavy  wall  is  considered.  So  far  as  the  author  is  aware,  this  is  the 
first  solution  obtained  for  a  multi- dimensional  problem  with  no  restric¬ 
tion  on  the  magnitude  of  the  absorption  coefficient.  The  results  sh  >w 
that,  in  general,  two  systems  of  waves  are  present  in  the  flow  field, 
a  modified  classical  wave  and  a  radiation- induced  wave.  The  occurence 
of  pressure  drag  at  subsonic  speeds  and  the  smoothing  of  the  transition 
from  subsonic  to  supersonic  speeds  reflect  the  nonequilibrium  character 
of  the  radiating  flow. 


II  FUNDAMENTALS  OF  RADIATION  TRANSPORT 

Most  of  the  existing  texts  on  thermal  radiation  vere  written 
by  astrophysicists  (for  example,  Chandrasekhar  (1950),  Koxrganoff  (1952), 
and  Sobolev  (1963)).  The  problems  with  which  they  deal  are  character¬ 
ized  by  an  absorbing,  emitting,  and  scattering  medium  that  is  stationary, 
time- independent,  and  infinit**  in  extent.  In  problems  of  radiative 
gaedynamics,  aerodynamic ists  must  often  be  concerned  with  time- dependent, 
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multi- dimensional  situations  in  the  presence  of  a  radiating  wall.  In 
this  chapter,  ve  will  review  briefly  the  fundamentals  of  radiation- 
transport  from  the  aerodynamic ist's  point  of  view.  A  useful  discussion 
on  various  aspects  of  radiation- transport  from  this  point  of  view  can 
also  be  found  in  the  works  of  Lighthill  (i960),  Vincent!  and  Baldwin 
( 1962 ) ,  and  Goulard  ( 1962 ) . 

2.1.  Photons . 

It  was  suggested  by  Planck  (1900)  that  the  emission  and 
absorption  of  radiant  energy  by  matter  does  not  take  place  continuously, 
but  in  finite  quanta  of  energy.  Einstein  (1905)  went  a  step  further  to 
suggest  that  not  only  is  the  radiant  energy  emitted  and  absorbed  in 
quanta  of  energy  but  that  it  travels  through  space  in  such  quanta  with 
the  speed  of  light.  The  quanta  of  energy  are  commonly  called  photons. 

It  is  understood  now  that  photons  have  the  dual  characteristics  of 
both  waves  and  particles.  They  show  wave  characteristics  with  regard 
to  propagation,  and  particle  characteristics  during  the  interaction 
between  radiation  and  matter.  When  the  wave  description  is  used, 
radiation  is  characterized  by  a  wave  length  X  and  a  frequency  v 
connected  by 


vX  =  c  ,  (2.1) 

where  c  is  the  speed  of  light.  When  the  particle  description  is 
used,  it  is  characterized  by  a  momentum  P  and  an  energy  E  connected 
by 
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(2.2) 


P 


where  is  the  unit  vector  in  the  direction  in  which  the  photon  is 
moving.  The  energy  E  in  Eq.  (2.2)  is  given  by  the  expression 


E  *  hv  , 


(2.3) 


where  h  is  Planck's  constant. 

2.2.  The  Distribution  Function  and  the  Radiation  Intensity. 

Since  photons  have  the  characteristics  of  particles,  we  may, 

as  in  the  kinetic  theory,  use  the  distribution  function  to  characterize 

the  radiation  field.  The  distribution  function  f  is  defined  as  the 

v 

number  of  photons  "belonging  to  dv  dft"  (i.e.,  photons  in  the  frequency 
interval  dv,  moving  within  the  element  of  solid  angle  dft  centered 
about  the  unit  vector  ft)  per  unit  volume  at  time  t. 

Astrophysicists  have  found  it  convenient  to  use  the  radiation 
intensity  to  characterize  the  radiation  field.  The  radiation  intensity 
is  defined  as  the  amount  of  radiant  energy  carried  by  photons  "belonging 
to  dv  dft"  that  crosses  a  unit  area  per  unit  time.  Consider  an  arbitrary 
area  element  dA  with  its  orientation  designated  by  a  unit  vector  n 
normal  to  dA.  (in  the  following,  we  shall  follow  the  convention  that 
the  orientation  of  the  surface  element  of  a  body  is  designated  by  the 
outward- drawn  unit  vector  normal  to  dA.)  During  a  time  interval  dt, 
those  photons  "belonging  to  dv  dft"  that  cross  dA  at  the  beginning  of 
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the  time  interval  have  travelled  a  distance  cdt.  The  trajectory  of 
these  photons  has  generated  a  cylinder  with  base  dA  and  slant  height 
cdt  in  the  direction  of  ft  Thus  the  number  of  such  photons  crossing 

dA  during  the  time  interval  dt  is 

f  c  ft  •  n  dA  dt  • 
v 

It  follows  that  the  amount  of  radiant  energy  carried  by  such  photons 
crossing  dA  during  time  interval  dt  i6 

hv  tv  c  ft  *  n  dA  dt  .  (2.5) 


3y  definition,  the  radiation  intensity  1^  is  obtained  by  dividing 
expression  (2.5)  by  the  time  interval  dt  and  the  area  element  per¬ 
pendicular  to  ft.  Ibus  we  have 

hv  f  c  ft  •  n  dAdt 

I  =  - 7 — 3 -  =  hv  c  f  .  (2.6) 

n  •  ft  dA  dt  v 

From  the  definitions  of  the  distribution  function  and  the  radiation 
intensity,  it  follows  that  these  quantities,  in  general,  are  functions 
of  frequency  v  (or  energy  E  *  hv),  position  ?,  direction  ft,  and 
time  t.  To  recognize  this  explicitly,  we  can  write 


fv  *  fv(r,  ft,  t)  , 
Iv  =  Iv(r,  ft,  t)  . 


(2.7) 
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Radiant  Heat-Flux  Vector,  and  Radiant  Pressure  Tensor. 

The  radiant  energy  density  is  defined  aB  the  radiant  energy 
per  unit  volume  per  unit  frequency  It  follovs  from  the  definition  of 
the  distribution  function  that  the  specific  radiant  energy  density  of 
photons  moving  in  the  direction  of  ft  is  hv  f  dft.  The  radiant  energy 

density  u  is  therefore  given  by 

v 

u  (r,  t)  =  /  hv  f  41  ,  (2.8) 

v  ft 

where  the  suffix  ft  signifies  that  the  integration  is  carried  out 
over  the  entire  range  of  solid  angle.  With  the  help  of  Eq.  (2.6),  Eq. 
(2.8)  can  be  expressed  in  terms  of  as 

ur  (t,  t)  -  £  /  Iy(r,  ft,  t)  dft  .  (2.9) 

v  ft 

Ihe  space- integrated  radiation  intensity  is  defined  by 

I  (f,  t)  =  /  Iy(r,  ft,  t)  dft  .  (2.10) 

v  ft 

It  follows  from  Eqs.  (2.9)  and  (2.10)  that  the  radiant  energy  density 
and  the  space- integrated  radiation  intensity  are  related  by 

IQ  (r,  t) 

ur  (r,  t)  =  — - -  .  (2.11) 
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We  are  now  In  a  position  to  find  the  radiant  heat-flux  vector 
and  the  radiant  pressure  tensor.  To  this  end,  let  us  suppose  that 
associated  with  each  photon,  there  is  some  property  g,  the  magnitude 
of  which  depends  on  c.  From  expression  (2.4),  it  follows  that 

gfycii  •  n  dA  dt  (2.12) 

represents  the  amount  of  this  property  transported  across  dA  during 
the  time  interval  dt  by  photons  "belonging  to  dv  d n".  Thus 

g  fv  c  a  •  n  dn  ,  (2.13) 

is  the  amount  of  g  per  unit  frequency,  carried  by  photons  moving 
within  the  solid  angle  d n  about  n,  per  unit  area  per  unit  time. 

The  net  flux  of  g  across  the  surface  element  is  obtained  by  integrating 
this  expression  over  all  directions.  Thus  we  have 

/  g  f  c  fi  •  n  dn  =  (n  •  /  g  f  c  n  dn)  *  (n  •  f)  ,  (2.14) 

n  n  v 

where  the  vector 


f  s/gfcn^n  ,  (?.15) 

n 

is  called  the  flux  vector  associated  with  the  property  g.  This  vector 
has  the  physical  significance  that  the  component  of  the  vector  in  the 
direction  n  is  the  flux  of  the  associated  physical  property  across  a 
surface  normal  to  n. 
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Consider  first  the  transfer  of  energy  flux.  If  we  let 
g  =  hv,  Eq.  (2.15)  becomes 

f«/hvcfnai-/Ifldn*  5  (r,  t),  (2.16) 

ci  Cl 

where  is  the  radiant  heat-flux  vector.  The  component  of  this 

vector  in  the  direction  parallel  to  the  coordinate  is 

q  (r,t)  =  n  •  q  =  /  n  •  0  I  (r,  ft,  t)  dft  =  /  l  I  (r,  Cl,  t)  d Cl, 

i  n  n 

(2.17) 

where  n.^  is  the  unit  vector  in  the  direction  x^  and  is  the 

direction  cosine  of  the  unit  vector  Cl  with  respect  to  x^  For  the 
coordinate  system  in  Fig.  1,  we  have  dft  =  sin  9  d0  d0,  and 
ft  =  /^n^  +  t  n^  +  l^y  where 

l l  ~  nl  ’  ^  =  sin  ®  6in  ^2  =  ^2  ^  *  COS  ^  ’ 

(2.18) 

=  n^  •  ft  =  sin  9  cos  0 

More  generally,  the  component  of  q^  in  the  direction  of  n  is  given 
by 

qv  (r»  t)  «  n  •  q  -  /  n  •  5  I  (r,  Cl,  -t)  d Cl.  (2.19) 

n  a 
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2.4.  Splitting  of  the  Space- Integrated  Itedlatlon  Intensity  and  the 
Radiant  Heat  Flux. 

It  is  sometimes  convenient  to  distinguish  quantities  associated 
with  photons  moving  in  certain  directions  by  superscripts  +  and  -. 
Imagine  a  transparent  plane  with  its  orientation  denoted  by  n.  The 
notations  f*  and  I*  are  respectively  the  distribution  function  and 
the  radiation  intensity  associated  with  those  photons  moving  in  the 
direction  of  the  unit  vector  where  fl  •  n  >  0.  The  space- integrated 
radiation  intensity  and  the  radiant  heat  flux  associated  with  these 
photons  are  given  respectively  by 
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(2.22) 


I*  (?,t)  =  /  I*(f,  <5,  t)  d(J, 

v  n-n>0 

and 


<1*  (r,t)  *  /  ft  .  n  I*(r,ft,t)  d ft,  (2.23) 

n  ft-n>0 

where  the  suffix  ft  •  n  >  0  signifies  that  the  integration  is  extended 
only  over  the  hemisphere  where  ft  •  n  >  0. 

Similarly,  if  f^  and  I~  are  respectively  the  distribution 
function  and  the  radiation  intensity  associated  with  photons  moving  in 
the  direction  of  the  unit  vector  ft  where  ft  •  n  <  0,  the  space- inte¬ 
grated  radiation  intensity  and  the  radiant  heat  flux  associated  with 
these  photons  are 


l‘  (r,  t)  *  /  f(r,  ft,  t)  dft  ,  (2.24) 

V  ft  *n<  0 

and 


q’  (r,  t)  *-  /  ft  ‘  n  l‘(?,  ft,  t)  dft,  (2.25) 

n  ft-n<0 

where  the  suffix  ft  •  n  <  0  signifies  that  the  integration  is  extended 

only  over  the  hemisphere  where  ft  •  n  <  0.  Note  that  I*  ,  q+  ,  I~ 

v  ‘’n  v 

and  q  given  respectively  by  Bqs.  (2.22)  through  (2.25)  are  all 
n 

positive  quantities. 
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It  follows  from  Eqs.  (2.10),  (2.19),  (2.22)  through  (2.25) 


that 


I0  (r,t)  =  Iq  (r,  t)  +  I*  (r,  t)  , 

V  V  V 

and 


Qv  (r,t)  =  q*  (r,  t)  -  q‘  (r,  t)  . 
n  n  n 


(2.26) 


(2.27) 


2.5.  Emission  Coefficient,  Absorption  Coefficient,  Scattering 
Coefficient,  and  the  Scattering  Function. 

Radiation  Interacts  with  matter  through  emission,  absorption, 
and  scattering.  The  absorption  of  photons  by  a  molecule  raises  its 
rotational,  vibrational,  or  electronic  energy  level.  Conversely,  by 
emitting  photons,  a  molecule  lowers  its  energy  level.  The  number  of 
photons  "belonging  to  dv  dfi"  emitted  in  the  volume  dV  during  time 
dt  is 


e 

(jjjJ)  dv  dC  dV  dt  ,  (2.28) 

where  e^  is  the  emission  coefficient  which  is  a  function  of  position, 
angular  direction,  frequency,  and  time. 

It  should  be  noted  that  photons  do  not  collide  with  each  other 
and  that  they  all  move  with  the  same  speed,  namely,  the  speed  of  light. 
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Therefore,  in  considering  the  collisions  between  photons  and  molecules, 
we  may  assume  that  the  molecules  are  at  rest  As  a  result  of  collisions, 
photons  may  be  absorbed  or  scattered.  The  decrease  in  photons  "belonging 
to  dv  d Cl"  in  volume  dV  during  time  dt,  as  a  result  of  absorption 
and  scattering  can  be  considered  to  be  proportional  to  the  path  length 
ds  and  the  number  of  such  photons  in  dV.  If  a  and  are  the 

absorption  and  scattering  coefficients  respectively,  the  quantities 

av  f  (r,  fl,  t)  dv  dfl  dV  ds  (2.29) 

and 

T|  f  (r,  fl,  t)  dv  dfl  dV  ds  (2.30) 

represent  the  decrease  in  such  photons  as  a  result  of  absorption  and 
scattering.  Both  ay  and  tj^  are  functions  of  frequency,  position, 
and  time.  The  reciprocal  of  the  absorption  coefficient  is  called  the 
mean  free  path  of  photons. 

It  follows  from  Eq,  (2.30)  that  the  number  of  photons 
"belonging  to  dv  dfl"  in  volume  dV  during  time  dt,  scattering  from 
the  direction  O'  to  other  directions  is 

tj y ( r , t )  f v(t,  fl*,  t)  dv  dfl  1  dV  ds  .  (2.31) 

Of  this  number,  a  certain  fraction  Ky(r,  fl'  -»  fl,  t)  dfl  will  emerge 
from  these  collisions  moving  in  directions  that  lie  in  the  element  of 
solid  angle  dfl  about  fl.  Here  Ky  is  called  the  scattering  function. 
Thus  the  expression 
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*v( r,  O'  -»fl,  t)  Tjv(r,  t)  f^(r,  fi'r  t)  dv  d O'  dV  ds  dfl,  (2.3?) 


represents  the  number  of  photons  initially  directed  along  fl*  which 
are  scattered  into  direction  about  ft  while  in  the  volume  dV.  The 
integral  of  this  expression,  namely, 

/  k  (r,  fl»  -»n,  t)  r)  (?,  t)  f  (r,  fl',  t)  dv  dtl*  dV  ds  dfl, 

(2.33) 

gives  the  total  number  of  such  contributions  from  all  directions.  The 
scattering  function  is  normalized  such  that 

/  k  (?,  6’  -»«,  t)  dnf  -  1  .  (2.3*0 

n' 


2.6.  Radiation-Transport  Equation. 

The  equation  of  radiation  transport  is  an  equation  of  con¬ 
servation  of  photons,  or  radiant  energy.  Let  us  consider  a  packet  of 
photons  "belonging  to  dv  dfi"  occupying  an  arbitrary  volume  element 
dV  at  time  t.  After  a  time  interval  dt,  the  number  of  photons  in 
the  packet  will  decrease  on  account  of  absorption  by  matter  and 
scattering  out  of  the  packet  as  the  result  of  collisions;  it  will  in¬ 
crease  on  account  of  emission  and  scattering  into  the  pqcket.  The  net 
rate  of  increase  of  photons  in  the  packet  is 
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(2.35) 


r 

r  Afv  -  1 

f  dv  dft  dV  « 

[  +  «  fl  •  grad  t  J 

dv  d a  dV. 


From  Eg.  (2.28),  the  rate  of  increase  of  such  photons  in  the  packet 
as  the  result  of  emission  in  dV  is 


e 

r-—  dv  dft  dV  .  (2.36) 


It  follows  from  Eq.  (2.33)  that  the  rate  of  increase  of  such  photons 
in  dV  as  a  result  of  collisions  which  scatter  photons  from  other 
directions  ft'  to  the  direction  fl  is 


t)  (r,-t)  ;  k  (f,  n'  -*a,  t)  f  (f,  n*, 

\j  V  V 

v  fl' 


t)  dv  dfl'  dv  ds  dfl/dt  . 

(2. 37) 


As  given  by  Eqs.  (2.29)  and  (2.30),  the  rate  of  decrease  of  such  photons 
as  a  result  of  absorption  and  scattering  is 


(av+Tjv)  fy(r,  fl,  t)  dv  dfl  dV  ds/dt  .  (2.38) 

From  conservation  requirements,  collecting  the  results  of  Eqs.  (2.35) 
through  Bq.  (2.30),  and  using  the  relation  ds  =  cdt,  we  have  finally 

«j-—  +  c  ft  •  grad  f^ 


(2.59) 


Eq.  (2.39)  is  the  radiation- transport  equation  in  terms  of  the  distribu¬ 
tion  function.  For  a  purely  scattering  medium,  this  equation  reduces 
to  a  form  similar  to  the  Boltzmann  equation  in  kinetic  theory.  The 
only  difference  lies  in  the  scattering  term.  In  kinetic  theory,  the 
integral  term  is  nonlinear  and  is  integrated  with  respect  to  velocity 
space,  whereas  in  the  radiation- transport  equation,  the  scattering 
term  is  linear  and  is  integrated  with  respect  to  solid  angle.  This  is 
due  to  the  fact  that  in  the  kinetic  theory  collision  is  between  like 
molecules  with  different  speeds,  whereas  in  radiation- transport  theory, 
collision  is  between  photons  with  constant  speed  and  molecules  at 
rest.  With  certain  changes  in  notation,  Eq.  (2.39)  Is  also  the  governing 
equation  for  neutron  transport  with  a  constant-cross-section  approximation. 

In  view  of  Eq.  (2.6),  EJq.  (2.39)  can  be  rewritten  in  terms 
of  radiation  intensity  as 

c  ST  +  a  •  gmd  Iv  (2. to) 

*  %  -  (<V\)  Iv  +  \  !  S'  -fi,  t)  Iv(?,  fi',  t)  «•. 

The  second  term  on  the  left-hand  side  of  Eq.  (2.40)  can  be  written  in 
a  number  of  different  forms.  For  example,  we  can  write 

bl  _  dl 

fi  *  grad  lv  =  *  div(fi  Iy)=  (2.41) 

V 

where  s  is  the  distance  in  the  direction  of  the  unit  vector  fi. 
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2.7.  Frequency- Integrated  Quantities . 


In  tne  foregoing  discussion,  we  have  used  the  subscript  v 
to  denote  monochromatic  quantities.  To  each  monochromatic  quantity 
Q^,  we  shall  now  associate  a  corresponding  frequency- integrated  quantity 
Q  defined  by 


Q  E  /  Q  dv  . 
0 


(2.42) 


Thus  the  expressions 


I  (r,t)  «  /  I  dv  =  /  /  I  dn  d v  -  /  l(r,n,t)  d n,  (2.43a) 

0  o  °v  o  n  v  n 


00  00 

u  (r, t)  *  /  u  d V  ■  -  /  /  I  dn  dv 


0  rv 


C  V 

on 


7  /  Kr,n,t)  dn, 
c  n 


(2.4>) 


q.(r,t)  2  /  q  dv 
0  i 


/  f  l  I  dn  dv 
i,  i  i  v 

o  n 


/  i±  i (?,  n,t)  an j 


n 


(2.43c) 


Pn(r,t)  =  /  P  dv  «  -  /  /  /  i  I  dn  dv 

1J  o  ij  o  n  J 


=  7  /  /./.  i(f,  n,  t)  dn,  (2.43d) 

c  n  " 


are  respectively  the  frequency- integrated  radiation  intensity,  radiant 
energy  density,  radiant  heat-flux  and  radiant  pressure. 


18 


2.8.  Complete  Thermodynamic  Equilibrium  and  Local  Thermodynamic 
Equilibrium. 


A  system  that  Is  simultaneously  In  mechanical,  thermal,  and 
chemical  equilibrium  Is  said  to  be  In  thermodynamic  equilibrium.  Thus 
for  a  system  to  be  In  thermodynamic  equilibrium,  the  temperature  and 
pressure  must  be  uniform  throughout  the  system  and  no  net  chemical 
change  can  take  place  within  the  system.  If  the  system  is  shielded 
from  external  radiation,  it  then  contains  uniform  density  of  photons 
moving  randomly  in  all  directions,  and  the  radiation  field  is  homogenous 
everywhere . 

If  we  denote  quantities  in  complete  thermodynamic  equilibrium 
by  superscript  *,  Bq.  (2.40)  gives,  in  view  of  Eq.  (2.34), 


e 

v 


a 


v 


* 


I 

v 


(2.44) 


It  is  shown  in  quantum  theory  that  if  the  system  is  in  complete  thermo¬ 
dynamic  equilibrium,  the  equilibrium  radiation  intensity  is  given  by 
the  Planck  function 


2hv' 


exp(|£)  -  1 


(2.45) 


where  K  is  the  Boltzmann  constant  and  T  is  the  temperature  of  the 
matter  as  well  as  of  the  radiation.  It  follows  from  Bqs.  (2.44)  and 
(2.45)  that  the  emission  and  absorption  coefficients  are  related  to 
the  Planck  function  as 
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(2.46a) 


ll  ,  2h^  _ 1 _ 

a  ~  2  .hVv  . 

V  C  expfe)  "  1 


In  astrophysics,  the  function  is  defined  as 


b  .II 


v  a. 


c  exp(^)  -  1 


(2.46b) 


For  a  radiation  field  in  equilibrium,  Eqs.  (2.9),  (2.10), 
(2.17)  and  (2.21)  yield,  with  the  aid  of  Eqs.  (2.18), 


4*  r  2hv?  1  “I 
c  5”  ,hv.  . 

L  C  CXp  KT  ‘  1  J 


(2.47a) 


I  —  4nB  , 
o  v 1 

v 


(2.47b) 


q  =  0  , 

vi 


(2.47c) 


pv  ,  *  ^  Bv  8ij  ' 


(2.47d) 


where  & ^ j  is  the  Kronecker  delta,  and  T  is  the  temperature  of  the 
gas.  It  follows  from  Eqs.  (2.47a)  and  (2.47d)  that  for  a  radiation 
field  in  equilibrium,  the  pressure  tensor  and  the  energy  density  are 
related  by 


*  i  # 

p  =  —  u  5,  . 
ViJ  5  rv  1J 


(2.48) 
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Eqs.  (2.22)  through  (2.25)  when  specialized  to  radiative  equilibrium 
become 


c  ■  C  *  2n  »v  *  (2-49,) 

v  v 

and 


q+#  =  q”#  =  irB  .  (2.49b) 

n  n 


The  frequency- integrated  Planch  function  is  given  by 


B 


/  B  dv 
0 


/ 

0  c 


2hv' 


dv 


/hvN 

eXP  KT 


-  1 


oT 

n 


(2.50a) 


where  o  *  (2ir  K4/l5h'c2)  is  the  Stefan- Boltzmann  constant.  If  we 
Integrate  Eqs.  (2.47)  through  (2.49)  over  the  complete  range  of  frequency, 
we  have,  in  view  of  Eq.  (2.50a) 


*  4oT4 
ur  =  "~c  * 

(2.50b) 

I*  =  4cTU  , 

0 

(2.50c) 

* 

q  =  0  , 

(2. 50d) 

n 

*  1  * 

4oT^ 

(2. 50e) 

PU  *  3  Ur  Blj  ' 

^  8iJ  ' 

4 

+*  -  #  oT 

1*1*  , 
0  0  If  ' 

(2.50f) 

+*  -*  4 

%  *  «n  ■  01  ' 

(2.50g) 
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In  tne  subsequent  discussion,  we  shall  define  a  black  surface  as  the 
surface  whose  I+,  I*,  and  q^  are  given  respectively  by  Eqs.  (2.50a), 
(2.50f),  and  (2.506). 

The  concept  of  local  thermo  dynamic  equilibrium  is  useful 
when  a  system  is  out  of  equilibrium.  We  will  speak  of  local  thermo¬ 
dynamic  equilibrium  if  at  every  point  of  the  matter  in  question  a 
definite  temperature  can  be  assigned.  It  should  be  noted  that  Eqs.  (2.46) 
holds  also  when  local  thermodynamic  equilibrium  is  assumed. 

With  the  assumption  of  local  thermodynamic  equilibrium,  Eq. 
(2.40)  becomes,  with  the  aid  of  Eq.  (2.46b), 

1 

—  +  ft  *  grad  Iv  (2.51) 

-  avBv-  •  (Vv'  Xv  +  \  /  -*5»  *)  Iv(f.  t)  •»’  • 


2.9  Grey-Gas  Approximation,  Asymptotic  Situations. 

For  most  engineering  problems,  the  effects  of  scattering  can 
be  neglected.  In  this  case,  the  integro- differential  equation  (2.51) 
reduces  to  the  purely  differential  equation 

.  dl  dl 

-  =  -a  (I  -  B  )  .  (2.52) 

c  at  ds  v'  v  v7  ' 

Integrating  Eq.  (2.52)  with  respect  to  the  entire  range  of  frequency, 
we  have 
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m  +  *  ‘  !0  “v(lv  •  V dv  •  (2-5J) 

The  grey-gas  approximation  consists  in  replacing  the  absorp¬ 
tion  coefficient  ay  by  a  constant  value  a  independent  of  frequency 
v.  With  this  simplification,  Bq.  (2.53)  becomes 

* ■°,<i‘b)  ’  (2-54) 

where  B  is  the  integrated  Planck  function  given  by  the  expression 
(2.50a).  Owing  to  the  very  large  mangitude  of  c,  the  unsteady  term 
in  the  radiation- transport  equation  is  always  much  smaller  than  other 
terms  for  aerodynamic  applications  and  can  therefore  be  neglected  even 
for  time- dependent  situations.  If  the  unsteady  term  is  neglected, 

Eq.  (2.5M  becomes 


g  =  -a(l-B)  .  (2.55) 

It  can  be  shown  that  the  general  solution  to  this  equation  is 

T/  x  -/  a(6)ds  .  ,  _/  tv  /  a(s)  ds  ,  _  -/  a(s)  ds 

I(s)  =  e  1  j  a(s')  B(s')  e'  ds'  +  C  e  1 

(2.56) 

where  C  is  to  be  determined  from  the  boundary  conditions.  If  we 
define  the  optical  thickness  as 
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(2.57) 


s 

x  =  /  a(s)  d s 
0 

and  let  the  boundary  value  of  I  at  x*  be  l(x*),  Eq.  (2.56)  becomes 

T 

I(x)  =  /  B(x')  exp[-(x-x' )]dx'  +  l( t*)  exp[-(x-x*)]  •  (2.58) 

x* 

Consider  now  the  asymptotic  expansion  of  Eq.  (2.58).  To 
this  end,  we  shall  speak  of  an  optically  thin  or  thick  gas  if  the 
optical  thickness  given  by  Eq.  (2.57)  is  much  greater  than  or  smaller 
than  unity.  Let  us  define  the  dimensionless  quantities  s  =  s //  and 
a  »  a/aQ,  where  2  is  a  characteristic  length  of  the  problem  and  a Q 
is  a  reference  value  of  a.  Equation  (2.57)  in  terms  of  these  dimen¬ 
sionless  quantities  becomes 

s  ^  ^ 

1  =  a  2  /  a  ds  .  (2.59) 

0 

Since  the  dimensionless  quantities  are  of  order  unity,  the  integral  in 
Eq  (2.59)  is  of  order  one.  It  follows  from  Eq.  (2.59)  that  the 
alternative  definition  of  optically  thin  or  thick  gas  depends  whether 
a  2  is  much  greater  than  or  smaller  than  unity. 

For  an  optically  thin  gas,  the  magnitude  of  the  quantities 
in  square  brackets  in  Eq.  (2.58)  is  much  less  Ahan  unity.  We  may 
therefore  expand  the  first  exponential  function  in  a  power  series  in 
terms  of  (x-x')  For  most  engineering  problems,  the  temperature  of 
the  wall  is  low  compared  with  that  of  the  gas  and  the  last  term  in 
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Eq,  (2.58)  can  therefore  be  neglected.  If  we  retain  only  the  first 
term  of  the  series,  Bq.  (2.56)  then  becomes 

T 

I  =  /  B(t')  dx'  .  (2.60) 

X* 


It  follows  that 


I  «  B 


(2.61) 


In  view  of  Eqs.  (2.41)  and  (2.50a),  the  radiation- transport  equation 
(2.550  in  this  case  can  then  be  written  as 


.  dl  aaT4 

5  -  —  • 


(2.62) 


That  is,  emission  is  a  first-order  effect  in  terms  of  the  absorption 
coefficient,  whereas  self -absorbing  is  of  second  order.  Thus  an  optically 
thin  gas  is  non-absorbing.  Equation  (2.62)  can  be  obtained  alternatively 
by  proceeding  directly  from  the  differential  equation  (2.55)  as  follows. 

We  first  expand  the  radiation  intensity  in  a  power  series  in  a,  namely, 

I  *  J]  a”  I%n^.  Substituting  this  series  into  Eq.  (2.55),  collecting 
n=l 

like  powers  of  a,  and  setting  the  coefficient  of  a  equal  to  zero, 
we  then  obtain  an  equation  which  is  identical  to  Eq.  (2.62). 

If  Bq.  (2.62)  is  integrated  over  the  entire  range  of  solid 
angle,  we  have,  with  the  help  of  Bq.  (2.43c) 
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(2.63) 


37J 


4aoT 


which  is  the  thin-gas  approximation. 

For  an  optically  thick  gas,  we  assume  that  the  function  B(t') 
in  Eq.  (2.58)  may  be  expanded  in  a  Taylor's  series  about  t'  =  t. 
Furthermore,  the  last  term  in  Eq.  (2.58)  is  taken  to  be  small  and  the 
lower  limit  of  integration  is  replaced  by  Thus  Eq.  (2.58)  becomes 


I(t) 


/  dx 


-(t-t* ) 


B(t)  +  (t'-t) 


dB(t) 

dr 


(t'-t)' 


d^B(x) 


+ 


dx 


=  B(x) 


<*B(t) 

dx 


d  B(x) 

+  -  *  -  •  •  • 
dx 


(2.64) 


Alternatively,  Eq.  (2.64)  can  be  obtained  from  the  differential  equation 
(2.55)  by  first  expanding  the  radiation  intensity  as 


I(s)  -  S  (i)"  I(n)(s)  .  (2.65s) 

n=0 

Substituting  this  series  into  Eq.  (2.55)#  collecting  like  powers  of  a, 
and  setting  the  coefficient  of  (l/a)n  equal  to  zero,  we  have 


(0) 


B 


.(1) 


•  (2) 


dl 


(0) 


ds 

(1) 


dl 


ds 


dB 

ds 

d2B 

IP 


(2.65b) 
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Substituting  Eqs.  (2.65b)  into  Eq.  (2.65a),  we  obtain  an  equation 
identical  to  Eq.  (2.64). 

If  we  retain  the  first  three  terms  in  Bq.  (2.6k),  and  sub¬ 
stitute  the  resulting  equation  into  Eq.  (2.43c),  we  have 


^(r) 


16oT^  dT 
3a  3x^ 


(2.66) 


The  radiant  heat-flux  given  by  Eq.  (2.66)  is  called  the  thick-gas  (or 
Rosseland)  approximation  It  is  in  a  form  similiar  to  that  for  heat 
conduction  with  a  variable  thermal  conductivity  given  by  k  =  (!6oT5/3a). 


2.10  Formal  Solution  of  the  One- Dimensional  Radiation-Transport  Equation. 

We  shall  now  restrict  ourselves  to  the  most  simple  class  of 
problems  in  which  the  flow  variables  depend  on  only  a  single  spatial 
coordinate,  say  Let  us  consider  a  semi- inf inite  extent  of  gas  on 

one  side  of  a  planar  solid  surface  perpendicular  to  the  coordinate  x^ 
and  located  at  x^  =  0.  If  we  define 


a(x2)  dx2  , 


(2.67a) 


it  follows  that  t  and  q  are  related  by 

t  l(»a>  - 


(2.67b) 


where  ^  *  l  *  cos  9. 
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If  we  assume  that  the  surface  is  black  and  has  a  temperature 


T  ,  the  value  of  I*  at  the  wall  is  given  by  Eq.  (2.50a)  as 


cT 

I  (0)  =  • 


(2.68) 


With  the  aid  of  relations  (2.67)  and  (2.68),  Eq.  (2.58)  can  then  be 
written  as 


I*(n,  M  )  - 1  ?  *V)  «p[-  ill  +  —2  exp[-  a  j  , 

(2.69b) 

which  is  the  radiation  intensity  directed  away  from  the  wall.  Similarly, 
if  nc  wall  exists  at  infinity,  the  radiation  intensity  directed  toward 
the  wall  as  given  by  Eq  (2.58)  is 

I'(ti,  M)  *  £/  t\v)  ^  j^-  ^  .  (2.69b) 


The  space- integrated  radiation  intensity  and  the  radiant  heat- 
flux  can  be  found  by  substituting  Eqs.  (2.69)  into  Eqs.  (2.4,5a)  and 
(2.43c).  This  leads  to 


T)  1  00 

I^(tj)  =2o/  T  (n  * )  E1(r)-rj')  dq'  +2  of  T  (tj'  Je^h'-h)  drj'42o  T^  ^(q), 

0  n 

(2.70) 


T)  1  »  ,  , 

q(n)  =  2o  /  T  (tj')  EL(h-ti')  dr|’  -  2o  /  T4(t)')  ^(rj’-qJdq'+aoT^  E,(t))  , 
0  r\  w  5 

(2.71) 


where  the  function  E  (z)  is  defined  by 

n'  '  * 

E  (z)  =  /  e’z/^  un‘2  du  .  (2.72) 

0 

This  function  has  the  properties 


E  (z)  .  -  E  (-z)  , 

n'  '  tl  * 

(2.7» 

En(z>  ‘  -  4  Wz>  > 

(2.73b) 

Enl0)  *  TTI’  lf  n  -  %  *•»  ••• 

(2.73c) 

En(-)  .  0  . 

(2.73d) 

Equations  (2.70)and  (2 . 71 )  are  the  solutions  obtained  by  Goulard  (i960) 
for  one- dimensional  problems  with  a  solid  black  boundary. 

III.  BASIC  GAS- DYNAMIC  EQUATIONS  AND  THE 
EXPONENTIAL  APPROXIMATION 

In  this  chapter,  we  shall  discuss  the  fundamental  equations 
of  radiation  gas-dynamics  under  the  following  assumptions: 

1.  The  gas  is  in  translational,  vibrational,  and  chemical  equilibrium. 
That  is,  the  effects  of  viscosity,  thermal  conductivity,  vibration, 
and  dissociation  are  neglected. 
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2.  The  gas  is  perfect. 

3.  The  gas  is  in  local  thermodynamic  equilibrium. 

4.  The  contributions  of  radiation  pressur'-  and  radiation  energy  are 
neglected. 

5.  The  effects  of  scattering  are  negligible. 

6.  The  absorption  coefficient  is  independent  of  frequency,  i.e.,  the 
grey-gas  approximation  can  be  applied. 

Assumptions  1  and  2  are  introduced  to  isolate  the  influence  of  radiation 
for  subsequent  discussion.  In  reality,  for  the  temperature  range  in 
which  the  effects  of  radiation  are  important,  the  gas  is  nonperfect  and 
other  nonequilibrium  processes  occur  as  well.  Assumptions  J>,  4,  and 
5  are  approximately  satisfied  for  most  engineering  applications  (see 
Lighthill  (i960))  Assumption  6  is  introduced  to  make  the  equations 
tractable. 

3-1  Basic  Equations. 

A  radiating  gas  can  be  considered  as  a  mixture  of  gas  molecules 
and  photons.  Since  photons  have  zero  rest  mass,  however,  the  continuity 
equation  for  the  radiating  gas  is  the  same  as  in  the  classical  case. 

Thus,  we  have 


where  p  is  the  mass  density  and  u^^  (i  =  1,  2,  3)  are  the  velocity 

components  The  derivative  D/Dt  =  d/dt  +  u  d/dx.  is  the  substantial 

J  J 

derivative,  where  the  repeated  dummy  subscripts  denote  the  summation 
convention. 
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If  radiation  pressure,  viscosity,  and  body  force  are  neglected 


the  momentum  equations  are 


(i  =  1,  2,  3)  ,  (3.2) 


where  p  is  the  pressure  of  the  gas. 

The  energy  equation  can  be  obtained  by  applying  the  first 
law  of  thermodynamics  to  a  fluid  element  of  unit  mass.  This  leads  to 


p  J*  -  52 

p  Dt  Dt 


0  , 


(3.3) 


where  q  is  the  component  of  the  radiation  heat-flux  and  h  is  the 
J 

specific  enthapy. 

For  a  gas  in  vibrational  and  chemical  equilibrium,  the  thermal 
and  caloric  equations  of  state  for  the  assumed  perfect  gas  are 


and 


p  *  pRT  , 


(3.M 


h 


(3.5) 


where  T  is  the  temperature  of  the  gas. 

The  radiation-transport  equation  for  a  non-scattering,  grey 
gas  in  local  thermodynamic  equilibrium  is  given  by  Eq.  (2.5*0  as 


dl 


,  _  oT 

-ad  -  — 


(3.6) 
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where  the  unsteady  term  has  been  neglected  for  the  reason  explained 
in  Section  2.9. 

The  radiation  field  is  coupled  with  the  flow  field  through 
the  integral 


t)  *  /  i(r,  n,  t)  J  an  .  (3-7) 

n 

Equations  (3*1)  to  (3-7)  are  the  governing  equations  for  a 
perfect,  inviscid,  nonconducting,  nonscattering,  grey-gas  in 
local  thermodynamic  equilibirum.  If  we  substitute  Eq.  (3  7)  into 
SQL  (3*3)>  it  is  apparent  that  the  governing  equations  are  of  integro- 
differential  form.  Exact  analytical  solutions  are  therefore  difficult, 
if  not  Impossible,  to  obtain. 

For  a  thin  or  thick  gas,  equations  (3-6)  and  (3  7)  are  replaced 
respectively  by  Eqs.  (2.63)  or  (2.66).  Either  one  of  these  equations 
together  with  the  gas-dynamic  equations  (3-1)  through  (3*5)  constitute 
a  determinate  set  of  purely  differential  equations.  Thus  in  these 
asymptotic  situations,  the  in tegro- differential  equations  reduce  to 
purely  differential  form. 

3.2.  One- Dimensional  Flow  Problems  and  the  Exponential  Approximation. 

For  one -dimensional  problems,  the  integral  terra  in  the 
differential  equations  can  be  removed  by  the  so  called  "exponential 
approximation".  With  this  approximation,  analytical  solutions  valid 
for  the  full  range  of  absorption  coefficient  are  possible  for  one¬ 
dimensional  problems  (see  Vincent!  and  Baldwin  (1962)).  Since  this 
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scheme  of  approximation  for  the  one- dimensional  problem  is  related  to 
the  spherical-harmonic  approximation  that  we  shall  discuss  in  great 
length  in  subsequent  chapters,  we  shall  discuss  the  exponential  approx¬ 
imation  briefly  here. 

For  one-dimensional  problems,  the  special  forms  of  the  gas- 
dynamic  equations  (3 - 1 )  through  (3-5)  are  obvious  and  we  need  not  write 
them  down  explicitly  here.  If  a  black  surface  with  temperature  is 

located  at  x0  =  0,  and  if  the  solution  of  Eq.  (3*6)  for  the  one-dimen¬ 
sional  situation  is  substituted  into  Eq.  (3-7)*  the  resulting  expression 
for  q  is  given  by  Bq.  (2.71)  as 

1  li  ®  k  u 

q(rj)  =  2o  /  T  ( TjO  EU(n-V)  <h'  -  2o  /  T  (q')  Ejq'-q)  dq  *  +OT  E  (q) 
0  q 

(3.8) 

The  exponential  approximation  consists  in  replacing  the  exponential 
integral  function  E^z)  in  Bq.  (3-8)  by  a  purely  exponential  -function 
of  the  form 


E^z)  =  -j  b2  e'bz  ,  (3.9) 

where  b  is  a  constant  and  cannot  be  uniquely  determined.  For  example, 
Vincenti  and  Baldwin  (1962)  take  b  *  1.562  whereas  Wick  (1964)  takes 
b  =  3/2.  As  we  shall  see  in  Section  4.3/  if  we  set  b  =  ,  the 

exponential  approximation  is  equivalent  to  the  later  spherical-harmonic 
approximation  when  that  approximation  is  specialized  to  one-dimensional 
problems  with  a  suitable  choice  of  boundary  conditions  at  the  wall. 
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With  the  exponential  approximation  given  by  Eq.  (3.9) ,  the 
exact  expression  (3.8)  can  now  be  approximated  by 

,  x  2b20  .  -b(n-rj * )  . 

q(l)  =  — r-  /  T  (V )  e  '  1  dtp 

?  0 

-  -£sp-  /  T4(q')  e"b(n'“Tl)  dq'  +  ^  o  T*  e"b\  (3-10) 

n 

If  Eq.  (3  10 )  is  differentiated  twice  and  the  undifferentiated  equation 
is  used  to  remove  the  integral,  it  can  be  shown  that  q  satisfies  the 
equation 


^-4  -  b2q  -  l6o  T5  ^  *  0.  (3-11) 

dr)  n 

This  equation  together  with  the  one- dimensional  form  of  Eqs.  (3*1) 
through  (3.5)  are  then  the  governing  equations  for  one- dimensional 
radiating  gas  flow  with  the  exponential  approximation.  Thus,  as  a 
result  of  this  approximation  the  exact  one- dimensional  lntegro- 
differential  equations  are  approximated  by  a  set  of  differential  equations 
valid  for  the  cumplete  range  of  absorption  coefficient.  Unfortunately, 
this  scheme  of  approximation  is  not  applicable  to  a  multi-dimensional 


radiating  gas 


IV  APPROXIMATE  RADIATION-TRANSPORT  EQUATIONS 


AND  BOUNDARY  CONDITIONS 


In  Section  2.6,  we  pointed  out  that  there  Is  a  mathematical 
anology  between  transport  of  neutrons  and  photons,  and  that  Eq.  (2.39) 
can  also  be  taken  as  governing  equation  for  transport  of  neutrons.  Much 
effort  has  been  devoted  by  nuclear  physicists  to  the  approximate  solution 
of  this  integro-differential  equation  in  connection  with  the  design  of 
nuclear  reactors.  The  most  useful  technique  for  this  purpose  is  the  so- 
called  "spherical- harmonic  method".  The  underlying  principle  of  this 
method  is  to  replace  the  exact  integro-differential  transport  equation 
by  an  infinite  set  of  coupled  equations  of  purely  differential  form. 

In  the  flow  of  a  radiating  gas,  if  the  scattering  term  is 
negligible  compared  with  other  terms,  the  radiation-transport  equation 
itself  is  already  in  purely  differential  form.  In  spite  of  this,  the 
complete  set  of  governing  equations  is  still  of  integro-differential 
form  because  the  radiation  heat  flux  in  the  energy  equation  and  the 
radiation  intensity  in  the  transport  equation  are  related  through  an 
integral  (see  Section  As  we  shall  see,  however  by  applying  the 

spherical- harmonic  method  to  the  radiation- transport  equation,  the  first 
approximation  of  the  resultant  equations  together  with  the  gasdynamic 
equations  constitute  a  determinate  set  of  purely  differential  equations 
valid  for  the  complete  range  of  absorption  coefficient. 
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4.1.  Approximate  Radiation-Transport  Equations. 

We  shall  now  obtain  the  approximate  radiation- transport 
equations  for  a  grey  gas  in  local  thermodynamic  equilibrium  by  means 
of  the  spherical-harmonic  method.  This  method  was  first  developed  by 
Chandrasekhar  (1944)  for  astrophysical  radiative-transport  problems 
with  plane  geometry.  It  was  later  studied  in  greater  detail  in  connec¬ 
tion  with  neutron- transport  problems  by  Mark  (1944),  (1945) »  Marshak 
(1945) ,  Eavison  (1957)>  and  Weinberg  and  Wigner  (1958)-  The  procedure 
of  the  method  is  also  similar  to  that  used  by  Grad  (1949)  in  kinetic 
theory. 

Integrating  Eq.  (2.51)  with  respect  to  the  whole  range  of 
frequency  and  assuming  a,  rj,  and  k  are  independent  of  frequency,  we 
have 


is  +n'-*rad  1 

*  221-  .  (a+rj)  I  +  tj  /  <c(r,  fl'  ->fi,  t)  I(r,  fi*,  t)  dft ' .  (4.1) 

O' 

Following  the  procedure  used  in  the  neutron-transport  theory,  we  expand 
the  radiation  intensity  in  a  series  of  spherical  harmonics  as  follows: 

oe  £ 

I(r,  fl,  t)  -  2  JJ  A/(r,  t)  Afi)  ,  (4.2) 

1=0  m =-£ 

where  the  A^(r,t)*s  are  functions  to  be  found  and  the  Y^(n)'s  are 
spherical  harmonics  (see  Appendix  A  for  a  brief  discussion  of  the  spherical 
harmonics).  We  assume  also  that  the  scattering  function  K  is  a  function 
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only  of  the  angle  between  the  directions  Q  and  fl ' .  We  call  this 
angle  Q  and  define  cos  6q  3  =  fl  •  fl ' .  Since  the  scattering 

function  is  a  function  of  the  single  variable  a  suitable  series 

representation  can  be  given  in  terms  of  the  Legendre  polynomial  P#(pQ). 
Thus  we  have 


k(t,  n '  ->n,  t)  =  k(u,  n*)  =  2  p#(cos  6  ) 

°  uo  1  1  ° 


(4.3) 


We  now  substitute  Eqs.  (4.2)  and  (4.3)  into  Eq  (4.1).  If  the  resultant 
equation  i6  multiplied  by  Y^(ft),  the  complex  conjugate  of  Y^(fi),  and 
integrated  over  the  entire  range  of  solid  angle,  Eq.  (4.1)  becomes, 
after  application  of  the  recurrence  and  orthogonality  relations  for 
the  spherical  harmonics  (see  Appendix  B  for  details). 


1  [  (i+m+2)  (i+m+1)]1/2  f  ^/+1  I 

c  3T  "  - 2TT7T31 - -  [  +  1  7)x“  J 


N  ,  .nl/2  r  & 

•m)  (i-m-l)]  '  l-l  . 

&(2£- ^ 


f] 


[(i-m+2)  (/-m+l)]1/2  T  ^/+1 

2(2 /+M  [ 


m-1 

rl 


( /  +m)  (/-Hn-l)]1/  2  T  ^Z- 1  ,  1  1 

- 2T27-1J  —  [  s; - J 


+  [(Z+tn+1)  (Z-m+1)]1'*  **£+ 1  +  [(/+m)  (i-m)]1  c  1 


2/-1 


^T 


,  .maoT4 
_TT_T  -  (a-Hj  A|  +  — 


60m  °0 Z  * 


(4.4) 
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where  5  and  bA,  are  the  Kronepker  delta.  Equation  (4.4)  represents 
Om  0/ 

an  infinite  set  of  differential  equations  for  an  infinite  number  of 
unknown  functions  a”.  This  set  is  completely  equivalent  to  the  integro- 
differential  equation  (4.1). 

We  can  obtain  an  N-th  approximation  to  the  foregoing  set 

(called  the  P  -approximation  in  neutron- transport  theory)  by  truncating 
N 

the  series  (4.2)  after  the  term  2=N.  It  is  known  from  neutron- transport 
theory  that  an  odd  numbered  approximation  is  more  accurate  than  the 
succeeding  even-numbered  approximation  and  that  the  first  approximation, 
the  P^-approximation  (often  called  the  diffusion  approximation),  is 
sufficiently  accurate  for  most  problems.  For  our  purpose,  therefore, 
we  limit  ourselves  for  the  present  to  the  first  approximation.  For  the 
first  approximation  (N=l),  Eq.  (4.4)  with  A1^  =  0  for  l  >  2  gives 


(4.5a) 


1,  m 


-1  :  c  ST-  + 


M 

Vs  l 


^0  1  4n  -1  ,  .  -1 

-  J  =  —  A1  -  (a+n)A1  , 


(4.5b) 
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1,  m  =  0  :  — 


1 

c  dt 


4n  .0 

T  i*1  Ai 


-  (a+tj)  A 


(4.5c) 


If  »  -  1  :  5 


*i 

3r 


vs 


[ 


*°o 


-  i 


*0]  k. 


K.l)  A. 


(a+rj)  Ax 

(4.5d) 


Thus  for  the  first  approximation,  if  T,  k^,  and  are  known,  we 

have  in  Eqs.  (4.5),  four  equations  for  the  four  unknowns  A^,  A~\  A^,  A^. 
The  quantities  and  are  obtained  as  follows  from  Eq.  (4.5)  by 

utilization  of  the  orthogonality  property  of  the  spherical  harmonics 
together  with  0q  (2.34): 


1 


/ 

SI' 


*(U0, 


«') 


YqV) 


an*  = 


l 


a 


K^0> 


Q  ' )  <*> 


1 

U7 


(4.6a) 


/  nQ  *(n0,  n*)  41  • 

-1  =  1^^,  *S(S')«(vfi')  *'  - r— — 

/  *(«„,  O')  «' 

O'  0 


» 


(4.6b) 


where  ^  is  the  average  cosine  of  scattering.  In  the  case  of  spherically 
symmetric  scattering,  p  *  0,  consequently,  *  0. 

Equations  (4.5)  can  be  rewritten  in  terms  of  the  space- 
integrated  radiation  intensity  and  the  radiation  heat-flux  as  follows: 

If  we  substitute  the  series  (4.2)  into  Eqs  (2  43a)  and  (2.43c)  and 
utilize  the  orthogonality  property  of  the  spherical  harmonics,  it  can 
be  shown  that 
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I0(r,  t)  =  4n  A °(r,  t)  ,  (4.7a) 

qx(f,  t)  «  -  i  -Jy*  ^A'i'(r,t)  +  A^(r,  t)J  ,  (4.7b) 

^(r,  t)  =  y  A°(r,  t)  ,  (4.7c) 

q^(r,  t)  =  ^A‘L(f,t)  -  A*(r,  t)j  .  (*-7<D 

In  viev  of  Eqs.  (4.6a)  and  (4.7),  Eq.  (4.5a)  can  be  rewritten  as 

c  sr  +  ^  *  •  a(Io  *  ‘t<rt’  >  -  (‘l-8a 


and  with  some  algebraic  manipulation,  Eqs.  (4.5b)  through  (4.5d)  can 
be  rewritten  as 

1  ^i  ^o 

c  3t“  +  "SF  =  "^a  +  1  qi>  (1=  (^-®b) 

Note  that  the  scattering  coefficient  enters  only  on  the  right-hand 
side  of  Eq.  (4.8b).  Thus  the  effect  of  scattering  is  equivalent  to 
a  modification  of  the  coefficient  on  the  right-hand  side  of  Eq.  (4.8b). 

If  the  unsteady  term  is  neglected  for  the  reason  explained 
in  Section  2.9,  the  approximate  radiation- transport  equations  (4.8)  for 
a  .ion- scattering,  grey  gas  in  local  thermodynamic  equilibrium  become 
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(4.9a) 


and 


dIo 


=  -  3aq1, 


(i  =  1,  2,  3)  . 


(4.9b) 


These  equations  can  be  rewritten  in  vector  form  as 


div  q  -  V  ■  q 


*Q(lo  -  4oT  ) 


(4.10a) 


and 


grad  Io  $  V  Iq  =  -3aq  . 


(4.10b) 


Equations  (4.9)  or  Eqs.  (4  10)  are  identical  to  the  four  approximate 
radiation- transport  equations  previously  obtained  by  the  author  (1964) 
by  means  of  a  moment  method.  If  Iq  is  eliminated  from  Eqs.  (4.9)  or 
Eqs.  (4.10),  the  resulting  equations  (except  for  sign  conventions  and 
notation)  are  identical  to  the  equations  cited  without  proof  by  Traugott 

(1963). 

With  the  aid  of  Eq  (4.7)  and  Eq  (A. 4)  in  Appendix  A,  the 
radiation  intensity  in  the  expression  (4.2)  for  the  first  approximation 
can  be  written 
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I(xi#  X2,  Xy  0,  e,  t) 


=  ^  Io(x1,x^  ,x^,t)  +  Sin  0  sin  e  +  ^  cos  6  +  cos  0  sin  6 

(4.11a) 

The  vector  form  of  this  equation,  in  view  of  Eq.  (2.18),  can  be  written 
as 


I(r,  fl,  t)  =  ^  [lQ(r,  t)  +  3q  •  oj  .  (4.11b) 

With  the  aid  of  Eq.  (4.10b),  this  equation  can  also  be  rewritten  as 

I(r,t5,t)  =  [l0(?,t)  -5  n  •  grad  ij  •  (4.11c) 

As  we  shall  see  in  Section  4.2,  either  one  of  Eqs.  (4.11)  is  required 
for  use  in  connection  with  the  boundary  conditions. 

It  is  of  interest  to  note  that  Eqs.  (4.9)  reduce  to  the 
usual  thin-  or  thick-gas  approximation  when  the  mean  free  path  of  the 
photons  (l/a)  is  much  larger  or  smaller  than  the  characteristic 
dimension  l  of  the  problem.  This  cat:  be  shown  by  converting  both 
the  dependent  and  independent  variables  in  Eqs.  (4.9)  to  dimensionless 
quantities  of  order  unity.  It  then  becomes  apparent  that  for  a/  «  1 
the  absorption  term  in  Eq.  (4.9a)  is  of  higher  order  than  the  rest  of 
the  terms.  Thus  for  al  «  1,  Eq.  (4.9a)  gives 
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On  the  other  hand,  if  al  »  1,  the  order- of -magnitude  estimate  shows 
that  the  term  on  the  left  in  Eq.  (4.9*)  is  negligible  compared  with 
other  terms.  This  leads  to 


I 


o 


(4.13) 


With  the  aid  of  this  result,  Eq.  (4.9b)  then  gives 

16oT3  dT 

qi  *  *  3a”  sj;  • 


(4.14) 


Equations  (4.12)  and  (4.14)  are  the  well-known  thin-  and  thick-gas 
(Rosseland)  approximations.  They  thus  agree  with  what  is  obtained  from 
asymptotic  expansion  of  the  solution  of  the  exact  radiation-transport 
equation  as  discussed  in  Section  2.9. 

For  problems  involving  only  one  spatial  dimension,  say  x^, 
the  corresponding  P^-approximation  for  the  radiation-transport  equation 
and  the  radiation  intensity  in  this  case  are  obtained  by  setting  q^ 
and  qz  equal  to  zero  in  Bqs.  (4.9)  and  (4.11).  These  results  can 
also  be  obtained  formally  by  returning  to  the  exact  radiation- transport 
equation  for  one- dimensional  problem  and  expanding  the  radiation  inten¬ 
sity  in  terms  of  Lagendre  polynomials.  Since  the  spherical-harmonic 
series  with  m  =  0  reduces  to  a  series  in  terms  of  Legendre  polynomials. 


we  can  carry  out  the  procedure  alternatively  by  putting  =  0  for 
m  ^  0  in  Eqs.  (4.2)  and  (4.4).  This  leads  to  the  series 


I(x,ii,t) 


2  Ai(x,t)  P/(p) 


(4.15) 


where  the  upper  index  0  on  the  right-hand  side  of  Eq.  (4.15)  and  the 
subscript  of  x  have  been  dropped.  The  infinite  set  of  radiation- 
transport  equations  for  problems  with  plane  geometry  become 

1  (i+1)  ^1+1  l  1 

c  dt  (2/+5)  dx  (2/-1)  dx 


A .  acT^ 

■  n^rrr  -  ^  h  +  *?-  5o*  •  (4-l6) 

For  the  first  approximation  (P^-approximation) ,  we  set 
A,.,  At,  ...  equal  to  zero  in  Eqs.  (4.15)  and  (4.l6).  The  resulting 
equations  can  then  be  put  in  terms  of  the  space- integrated  intensity 
and  the  radiation  heat  flux,  which  are  given  by  Bqs.  (4.7a)  and  (4.7c) 

as 


IQ(x,  t)  =  4n  Ac(x,  t)  ,  (4.17a) 

and 

q(x,  t)  «  4r  A^(x,  t)>  (4.17b) 

where  the  subscript  of  q  in  Eq.  (4.17b)  has  been  dropped. 
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In  view  of  Bis,  (4.6)  and  (4.17),  the  P^-approximation  of 
Bis.  (4.15)  and  (4.16)  are 

I(x,  u,  t)  «  ^  [IQ(x,t)  ^  3u  q(*,t)],  (4.18) 


and 


i  x-. 

-  ^  -  “^(Iq  ”  4oT  )  ,  (4.19a  t 

,  ^  hi 

—  ^  *  -3(a+q-kiTi)  q  .  4.19b/ 


If  the  unsteady  term  and  the  effect  of  scattering  are  neglected.  Bis. 
(4.19)  become 


-a(lQ  -  401*^) 


(4.20a) 


=  -3oq  .  (4.20b^ 

Biuations  (4.20)  are  the  one -dimensional  case  of  Bis  (4.9).  Except 
for  sign  convention  and  notation,  these  equations  are  identical  to  those 
obtained  by  Traugott  (1963)  for  the  first  approximation  by  means  of  a 
moment  method.  The  PT-approximation  of  Bq.  (4.16),  neglecting  the 
effect  of  scattering,  is  essentially  the  second  approximation  obtained 
by  Traugott  (1963). 
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4.2.  Approximate  Boundary  Conditions. 

If  a  surface  Is  present  In  the  flow  field,  the  surface  condition 
and  the  distribution  of  wall  temperature  will  effect  the  radiation  and 
flow  fields  Suppose  a  black  surface  with  outward-drawn  unit  normal 
vector  n  is  described  by  the  equation  F(r)  =  0.  The  radiation 
intensity  at  the  wall  is  therefore  given  by 

oTj(?,t)  _  _ 

I  (r,  fl,  t)  =  — -  for  ft  •  n  >  0  ,  (4.21) 

n 

where  r  is  on  the  prescribed  surface  and  is  the  local  wall 

temperature  When  the  temperature  field  is  known,  the  exact  solution 
for  I  is  uniquely  determined  by  Eq.  (4.1)  subject  to  boundary  condition 
(4.21)  If  the  spherical-harmonic  method  is  applied,  the  integro- 
differential  equation  (4.1)  is  now  replaced  by  an  infinite  set  of  purely 
differential  equations  in  terms  of  A™  as  given  by  Eq.  (4.4).  The 
corresponding  infinite  set  of  boundary  conditions  in  terms  of  A®  can 
be  obtained  by  substituting  Bq.  (4.2)  into  Eq.  (4.21). 

The  infinite  number  of  conditions  as  given  in  Eq.  (4.21)  can¬ 
not  be  exactly  satisfied  in  an  approximation  of  finite  order.  Thus  for 
the  P  -  approximation,  another  set  of  approximate  conditions  consistent 
with  the  spherical-harmonic  method  will  have  to  be  found.  This  approx¬ 
imate  set  of  conditions  is  not  uniquely  defined.  The  most  frequently 
applied  boundary  conditions  are  those  which  were  proposed  by  Mark  (1944) 
and  Marshak  (1947).  Both  of  these  boundary  conditions  were  originally 
proposed  for  neutron-transport  problems  with  plane  geometry.  We  shall 
now  consider  the  application  of  Mark's  and  Marshak's  boundary  conditions 
to  radiation-transport  problems  with  plane  geometry. 
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In  the  special  case  of  a  one- dimensional  problem  with  a  black 


wall  located  at  *  0,  the  boundary  condition  (4.21)  can  be  written 


I  (0,  n,  t) 


oT4(t) 

w 


for  u  >  0 


(4.22) 


where  u  ■  cos  0.  Mark  proposed  that  for  the  p^-approximation  this 
condition  is  satisfied  for  certain  discrete  values  of  u  that  are 
taken  to  be  the  roots  of  the  Legendre  polynomial  *  0.  Mark's 

boundary  condition  in  the  P^-approximation  is  thus  written 

+  oT'(t) 

I  (0,  t)  *  — 2 -  for  ^  >  0  and  given  by  PN+1(^i)*°» 

(4.23) 

where  I*  on  the  left-hand  side  is  to  be  represented  by  the  truncated 
series  (4  15)-  We  note  that  Eq.  (4.23)  tends  to  the  exact  boundary 
condition  (4.22)  as  N  approaches  infinity. 

Marshak  suggested  alternatively  that  for  the  PN -approximation, 
the  boundary  condition  (4  22)  be  approximated  by 

1  +  1  oT 1 

/  I  (0,  H,  t)  P21_1(m)  (t)  P21-1(n)  du 

(i  *  1,  2,  ...  ,  |  (N+l) ) ,  (4.24) 

with  I  represented  by  the  truncated  series  (4.13)  • 


We  now  discuss  the  P^-approxlmation  of  Mark's  and  Marshak '  s 
boundary  conditions  For  the  P^-approximation,  the  radiation  intensity 
for  problems  with  plane  geometry  is  given  by  Eq ,  (4.18).  Substituting 
this  truncated  series  into  the  left-hand  side  of  0^.  (4.23)  and  noting 
that  the  positive  root  of  Pg^)  =  ^i  "  ^  ~  ^  is  V  V3  >  we  have 

if  _  “I 

]“  I  Io(0,  t)  +  V5  q(0,t)l  =  — 1 -  ,  (4.25) 

which  is  Mark's  boundary  condition  for  the  P^-approximation. 

For  the  P^-approximation,  the  left-hand  side  of  Eq.  (4.24) 
represents  the  heat  flux  emitted  from  the  wall.  If  we  substitute  Eq. 
(4.18)  into  Eq.  (4.24)  and  evaluate  the  integrals,  we  have 

l  [Io(0,  t)  +  2q(0,  t)]  *  °Ty ( t ) #  (4.26) 

which  is  Marshal's  boundary  condition  for  the  P^-approximation. 

Equations  (4.25)  and  (4.26)  can  be  both  written  in  the  form 

Io(0,  t)  +  mqfO,  t)  =  4oT^(t),  (4.27) 

where  m  =  V3  for  Mark's  boundary  condition  and  m  =  2  for  Marshak's 
boundary  condition. 

It  is  evident  from  the  foregoing  discussion  that  in  the  first 
approximation,  Mark's  boundary  condition  specifies  the  radiation  intensity 
to  be  satisfied  at  a  specific  value  of  p,  i.e.,  at  p  =  1/7/3  whereas 
Marshak's  boundary  condition  specifies  the  heat  flux  emitted  from  the  wall. 
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It  has  been  found  in  neutron-transport  theoiy  that  in  a  given  approximation 
the  difference  in  results  obtained  by  applying  Mark's  or  Marshak's 
boundary  conditions  is  insignificantly  small. 

*fhile  the  extension  of  Mark's  boundary  condition  to  multi¬ 
dimensional  problems  is  not  apparent,  the  generalization  of  Marshak's 
boundary  condition  leads  to  further  complications  and  ambiguities.  For 
a  general  treatement  of  this  matter,  the  reader  is  referred  to  Davison 
(1957 >  p.  167).  For  our  purpose,  it  is  sufficient  to  quote  without  proof 
the  P^-approximation  of  the  generalized  Marshak's  boundary  condition. 

To  be  specific,  consider  a  planar  wall  located  at  x^  *  0. 

The  radiation  intensity  at  the  wall  i6  then  given  by 

oT4 

I  O  s  (4.28) 


The  P^-approximation  of  the  Marshak-type  boundary  condition  becomes 

n/2  0=2  it  n  n/2  0*2n  oT^  _ 

/  /  I  Y°(0)  dO  =  /  /  Y?(0)  dO  .  (4.29) 

0=0  0=0  0=0  0=0 


In  view  of  Eqs.  (2.23)  and  (A-4a)  in  Appendix  A,  the  left-hand  side  of  Eq. 
(4.29)  represents  the  heat  flux  emitted  from  the  wall.  Substituting 
Bq.  (4.11a)  into  Eq.  (4.29)  and  evaluating  the  integrals,  we  have 


^  [Iq^,  0,  Xy  t)  +  2<i2(xl,  0,  Xy  t )  ]  «  °T^(x^»xj, t) . 

(4.30) 

This  is  the  P^-approximation  of  the  Marshak's  boundary  condition  for 
multi-dimensional  problems.  Comparing  Eqs.  (4.26)  and  (4.30),  we 
note  that  these  equations  are  of  the  same  form  and  that  Eq.  (4.26) 


is  a  special  case  of  Eq.  (4.30).  Equation  (4.30)  provides  a 
suitable  radiation  boundary  condition  for  the  planar  wall  in  a  three* 
dimensional  flow  field. 

4.3*  Hie  Relation  of  the  Exponential  Approximation  and  the  Spherical- 
Harmonic  Approximation. 

It  is  shown  in  this  section  that  the  exponential  approximation 
discussed  in  Section  3*2  is  equivalent  to  the  P^-approximation  when 
that  approximation  is  specialized  to  one-dimensional  problems  and  used 
with  Mark's  boundary  condition.  Ve  are  concerned  here  specifically  with 
the  one- dimensional  problem  posed  In  the  previous  section  with  boundary 
condition  given  by  Eq.  (4.22).  To  begin,  we  note  that  if  Iq  is 
eliminated  from  Eqs.  (4.20)  the  resulting  equation  Is  identical  to 
Eq.  (3-11)  with  b  =  V3  .  This,  however,  is  not  sufficient  to  prove 
that  the  two  schemes  of  approximation  are  equivalent,  since  the  informa¬ 
tion  on  boundary  conditions  is  lost  in  the  exponential  approximation 
when  the  integral  equation  is  converted  into  differential  equation. 

In  order  to  prove  the  equivalence,  we  first  eliminate  q 
from  Eqs.  (4.20)  to  obtain 

d2 10  „ 

- ^  -  3I0  ♦  12 or  =  0  .  (4.31) 

drj 

The  general  solution  of  this  equation  can  be  shown  to  be 

iQ(n)  =  2V3 of1  t\v)  (n,‘T,)  dn?  (4.32) 

0 

+  2Y3o  /  T*(V)  dT).  +  e*  V3n  4  f 

n 
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where  C^(t)  and  C(t)  are  to  be  determined  from  the  boundary  conditions. 
If  Iq  is  finite  at  infinity,  C2(t)  must  vanish.  Substitution  of 
Eq.  (4.52)  with  C2(t)  «  0  into  Eq.  (4.20b)  yields 

q(q)  =  2 0  /  T4(n’)  ^  dij' 

0 

-  2o  /  tNh')  +_i_  e*V3i  _  (4.33) 

n  V? 

Imposing  Mark's  boundary  condition  (4. 25)  and  utilizing  Bqs.  (4.52)  and 
(4.33),  ve  find 


CL(t)  -  2oT4  .  (4.34) 

Equations  (4.32)  and  (4.33)  together  with  expression  (4.34)  provide  a 
solution  for  I  and  q  as  obtained  from  the  P^-approximation 
specillzed  to  the  one- dimens  Iona  I  problem  and  utilizing  Mark's  boundary 
condition.  If  we  compare  these  results  with  the  exact  solution  given 
in  Eqs.  (2.70)  and  (2.71),  it  is  evident  that  the  Pj-approximation 
with  Mark '8  boundary  condition  is  equivalent  to  putting  b  «  V3  and 
hence 


EgU)  -  ,  (4.35) 

in  the  exact  solution  for  the  one-dimensional  problem.  It  should  be 
emphazled  that  our  conclusion  is  rather  general  and  is  Independent  of 
the  temperature  distribution  in  the  flow  field. 
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V.  LINEARIZED  THEORY 


It  was  shown  In  Section  4.1  that  to  a  first  approximation 
the  exact  multi- dimensional  transport  equation  can  be  replaced  by  four 
differential  equations  in  terms  of  the  unknowns  I  and  q^.  The 
approximate  radiation-transport  equations  given  by  Eqs.  (4.9)  together 
with  the  gas  dynamic  equations  (3*1)  through  (3*5)  constitute  a  deter¬ 
minate  set  of  purely  differential  equations  valid  for  the  complete 
range  of  absorption  coefficient.  Although  the  governing  equations  are 
in  differential  form,  they  are  nonlinear,  and  analytical  solution  is 
therefore  still  difficult.  In  this  chapter,  we  shall  examine  the 
corresponding  linearized  equations. 

5.1.  Acoustic  Equations. 

We  are  concerned  here  with  the  problem  of  weak  disturbances 
propagating  into  a  uniform  medium  that  is  in  radiation  equilibrium.  If 
the  coordinate  system  is  fixed  in  the  undisturbed  fluid,  the  dependent 
variables  are  given  by  u.  =u',p  =  p  +p,,p  =  p_+p,,a  =  a  -K*', 

I  =  I  +  1^,  qj^  =  q^  etc.,  where  the  disturbance  quantities  are 
denoted  by  primes  and  the  subscript  *  denotes  the  equilibrium  reference 
condition.  Substituting  these  variables  into  the  governing  equations 
(3.1)  through  (3»5)  and  (4.9)#  and  proceeding  with  the  linearization  in 
the  usual  fashion,  we  obtain 
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(5.1) 


q. 


sr 


0 


(i  -  1,  2,  3),  (5.2) 


ah' 

p-  sr 


(5.3) 


T' 


i 


(5.4) 


h' 


(5.5) 


ajl’  -  16oT^  T')  , 


(5.6a) 


3ol  q< 


(i  =  l,  2,  3)  •  (5.6b) 


We  now  Introduce  a  perturbation  velocity  potential  0  such 


that 


and 


P' 


(5.7) 
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These  expressions  satisfy  the  momentum  equations  (5*2).  The  remaining 
equations  (5.1)  through  (5*6)  can  then  be  reduced  to  a  single  equation 
in  terras  of  0.  To  this  end,  we  first  differentiate  Eq.  (5*^)  with 
respect  to  t  and  make  use  of  Eqs.  (5*1)  and  (5*7)  to  obtain 


dT'  1  f  d20  _2  d20  1 

=  ‘  *  L  S?  T 


(5.8) 


wheie  a,p  =  "YRT^  is  the  isothermal  speed  of  sound.  Eliminating  h* 
from  Eqs.  (5*3)#  (5-5)  and  making  use  of  Eqs.  (5.1)>  (5. 7)  and  (5.6a), 

we  have 


Poo 

7-1 


.2  a H  I 

sST5ijJ 


=  -  ajl°  -  16otV), 


(5.9) 


where  ag  s  Y?RT^  1®  the  isentropic  speed  of  sound.  If  q^  is 
eliminated  from  Eqs.  (5-6),  we  obtain 

V  3cr  I’  +  48a  otV  =  0  .  (5.10) 

J  00  O  00  00 


a  (K 


Finally,  Eqs.  (5.8),  (5*9)  and  (5*10)  are  combined  to  give  the  following 
single  equation  in  terms  of  0: 


i67«sa„  .S&S 

dt  dx.  dXj  +  Bo  dXj  dXj  "  ^  * 


(5.11) 


where  W„ 


7*0.  «S 

(r-i)  ot^ 


a20 


and 


Bo  = 


The  dimensionless  quantity 


Bo  is  called  the 


Boltzomnn  number.* 


Equation  (5-11)  is  a  fifth* order  linear  partial  differential 
equation  in  0.  With  appropriate  boundary  conditions,  0  can  be  deter¬ 
mined.  The  disturbed  velocity  and  pressure  fields  can  then  be  found 
from  Eqs.  ( 5 -  T )  •  lbs  disturbed  radiation  field  in  terms  of  0  can  be 
found  as  follows:  Eliminating  1^  from  (5.6b)  and  (5*9)  and  differen¬ 
tiating  the  resulting  equation  with  respect  to  t,  we  have,  in  view 
of  Bq.  (5.8), 


P-  AS  .  16°^ 

Differentiating  Eq.  (5.6a)  with  respect  to  t,  we  have,  in  view  of 
Eqs.  (5.8)  and  (5-12), 

dl'  ,  P_  16oT* 

W2  -  16 <fl*  - , -  ■;  - *—  .  L  .  (5.15») 

With  the  aid  of  Bq.  (5.11),  this  equation  can  also  be  written 

dIo  „jap  dws  l6pT-5  „  p- 

3t~  *  16oT-  3T  *  a.(r-i)  3T  1  '  ~jT  WT  *  ajr-l)  ST  • 

(5.15b) 

As  we  shall  see,  Eqs.  (5*12)  and  (5.13)  are  needed  in  formulating  the 
boundary  condition  at  the  wall. 

* 

By  using  the  exponential  approximation  applicable  in  the  one-dimensional 
case.  Lick  (1964)  obtained  au  equation  that  corresponds  to  Bq.  (5.11) 
specialized  to  unsteady  one- dimensional  flow. 


Equation  (5.11)  plays  the  same  role  here  as  does  the  wave 
equation  in  classical  acoustic  theory.  Since  Bq.  (5*11)  is  of  higher 
order  than  the  classical  equation,  additional  boundary  conditions  are 
necessary.  If  a  radiating  wall  exists  in  the  flow  field,  the  surface 
condition  as  well  as  the  distribution  of  wail  temperature  will  affect 
the  flow  field. 

To  simplify  the  analysis,  we  assume  that  the  wall  is  black. 
Since  viscosity  and  heat  conductivity  are  neglected,  the  temperature 
of  the  wall  is  not  necessarily  equal  to  the  temperature  of  the  gas 
immediately  adjacent  to  the  wall.  Consider  a  planar  black  wall  with 
local  temperature  Tw  located  at  xg  =  0.  The  relation  governing  the 
temperature  Jump  at  the  wall  can  be  obtained  by  utilizing  Eq.  (4.50) 
as  follows:  Linearization  of  Eq.  (4.50)  leads  first  to  the  form 

16oT^  T^(xl,  Xy  t)  *  0,  Xy  t)  +  2q^(x1#  0,  Xy  t).  (5.14) 

Differentiating  this  relation  with  respect  to  t  and  substituting 
Eqs.  (5-12)  and  (5*  15b)  into  the  resulting  equation,  we  obtaiu  finally 

r  K  i 

[  3 r  (V  xy  •  3T  (X1>  °'  *J'  *>J 


Bo  3  I"  _2_  dWS 

1^5^  3t  [  3a.  3^ 


(x1#0,x5,t) 


(xL,0,x5,t) 


(5.15a) 
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This  equation  shows  clearly  the  temperature  Jump  at  the  wall.  For  an 
opaque  gas  (aM -»«•),  Eq.  (5*15aj  becomes 


(x^#x^#t)  =  (x^,  0,  Xji  t)  .  (5*15b) 


Thus  for  an  opaque  gas  no  temperature  Jump  exists  at  the  wall. 

In  view  of  Bq.  (5.8),  Eq.  (5.15a)  can  also  be  written  as 


K 

3T  (*i»  xy  *> 


bo  a  r 

'  2  „  I  A  1  [ 

‘  2  , 

.  l 

"  l&V'Big  3t  [ 

3oT  s:  -  wsJ  +  B 

■  "  2  J(X  ,0,x  ,t)  1 

TJ 

(5.15c) 

This  provides  the  radiative  boundary  condition  at 

the  wall. 

5.2.  Galilean  Transformation;  Steady  Flow. 

For  problems  Involving  a  body  moving  with  constant  speed  and 
attitude  In  a  medium  at  rest.  It  is  convenient  to  Introduce  a  coordinate 
system  fixed  in  the  body.  In  such  a  coordinate  system,  the  problem 
transforms  to  one  of  steady  flow  about  a  stationary  body.  If  the  boty 
Is  slender,  a  single  equation  In  terms  of  0  can  be  obtained  by  returning 
to  the  governing  equations  (3.1)  through  (3*5)  and  (4.9)  and  proceeding 
with  the  linearization  accordingly.  The  desired  results  can  also  be 
obtained  from  the  preceding  section  ,  however,  by  a  Galilean  trans¬ 
formation.  If  the  body  is  moving  with  constant  speed  Uw  in  the 


positive  x^-direction,  Bq.  (5.11)  after  such  a  transformation  becomes 


1  J 


16a.  a  Lt  2  dLg 
~  5*  Ax  '  **m  33^  -  0  ' 


(5.16) 


■  2j  ^2j  ,2. 


where  Lg  *  (1-m|  )  +  0  +  0 , E  (1-m£  )  ||  +  M  +  M  ' 

00  ox^  dXp  dx^  »  dx^  dXg  ox^ 


^Rp.  U«e 

and  Bo  - - r  .  Here  M„  *  U  /a_  is  the  isentropic  Mach  number, 

(7-l)<  S«  *  S« 

end  M,  ■  U  J*  -  tM  is  the  isothermal  Mach  number.  With  the 

^oo  no  (jo 

transformation,  Eqs.  (5*7),  (5.8),  (5.12)  and  (5.15b)  become 


(5.17a) 


p'  =  -  p  U 

r  ^ao  t 


«  dx 


(5.17b) 


(5.17c) 


ST 


e.  uf  &% 

Jaf(r-i)  m!  dxi 


16oT5  U  dLm 

00  00  T 

*V»4 


(5.17d) 
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aio  j  dr  .  p-  u-  dLs 
^  =  16<  3^  *  qJr"irM7  ^ 


16oT5  U 

m  m 


L_  ♦ 
T 


P.  aLB 

a.(r-i)»4 


(5-ITe) 


After  the  Gellleen  trensforaetlon,  Etjs.  (5.15»)  end  (5*15c)  become 


[  3^  (xi*  xj>  *  (xi>  °>  xj)] 


Bo  U 


16a 


/HMg 


*[■ 


It  il' 


_2_  dLsl 

3a. 


(x1#0,x3) 


2Um  dLy 

— ~~2~  sr  (*i*°^x3)  > 


(5.1Q») 


and 


dT' 

srf  <v  V 


Bo  U 
00 

i&r^f 


_2_  dLSl 
3a.  3^J 

2  (x^OjX^) 


+ 


[lt  -  -k  i] 


(xL,0,x5) 


(5.18b) 


Equation  (5-16)  plays  the  same  role  in  radiating  gas  flow  as 
does  the  Prandtl-Glauert  equation  in  classical  compressible  flow  theory. 

As  we  recall,  for  supersonic  flow  in  classical  theory,  the  structure  of 
the  Prandtl-Glauert  equation  in  two  or  three  space  corrdinates  is  similar  to 
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that  of  the  acoustic  equation  in  one  less  space  coordinate.  The  streamwise 
coordinate  in  the  Prandtl-Glauert  equation  is  then  "time-like."  From 
the  works  of  Vincenti  (1959)>  Moore  and  Gibson  (i960),  and  Chu  (1957)# 
such  analogy  also  exists  in  flow  with  chemical  and  vibrational  non¬ 
equilibrium.  The  analogy  does  not  hold,  however,  in  the  flow  of  a 
radiating  gas.  This  can  be  seen  by  writing  Eq.  (5- 11)  explicitly  for 
the  two  space  coordinates  x^  and  x^  and  comparing  the  result  with 
Eq.  (5.16)  for  the  three  space  coordinates  x^,  and  Xy  This  leads 
to  a  comparison  of 
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The  structure  of  these  equations  is  not  similar  owing  to  the  additional 
terms  with  d  /bxl  in  Eq.  (5.19b).  The  reason  that  the  analogy  does 
not  exist  in  the  radiating  gas  flow  is  that  a  signal  in  the  radiation 
field  can  propagate  with  the  speed  of  light.  Consequently,  disturbances 
will  propagate  in  the  negative  direction  in  steady  flow  even  against  a 
supersonic  stream.  On  the  other  hand,  in  unsteady  flow  no  signal  can 
propagate  in  the  negative  direction  in  time. 
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VI.  ILLUSTRATIVE  EXAMPLES 


To  Illustrate  the  application  of  the  present  formulation  of 
radiating  gas  flow,  ve  shall  consider  several  simple  examples.  In 
particular,  we  will  reconsider  the  problem  of  one- dimensional  wave 
propagation,  which  has  been  solved  previously  on  the  basis  of  the 
exponential  approximation.  To  provide  a  simple  example  of  a  multi¬ 
dimensional  radiating  gas  flow  with  no  restriction  on  absorption  co¬ 
efficient,  we  shall  also  consider  the  steady  flow  over  a  wavy  wall. 

6.1.  Propagation  of  Linearized  One- Dimensional  Waves. 


In  this  section,  we  shall  consider  small  disturbances  generated 

by  a  planar  wall  perpendicular  to  the  x-axis  and  translating  normal  to 

itself. 

Specializing  Eq.  (5.11)  to  the  present  problem,  we  have 

l6r*S  a-  .2*8  „ 

bo  yr-30-^0  • 

(6.1) 

The  boundary  conditions  imposed  by  the  motion  and  temperature 

of  the  wall  are 

variation 

u'(0,  t)  =  ^  (0,  t)  »  given  function  of  t. 

(6.2a) 

T^(t)  =  given  function  of  t. 

(6.2b) 

6l 

The  boundary  condition  at  infinity  is 


0(«,  t)  *  finite  quantity.  (6.3) 


The  boundary  condition  (6.2b)  can  be  expressed  in  terms  of 
0  as  follows:  Since  the  approximate  radiation  boundary  condition  (4.26) 
for  one-dimensional  problems  is  a  special  case  of  Eq.  (4.30),  we  can 
specialize  Eqs.  (5.15a)  and  (5*  15c)  for  one- dimensional  problems  to 
obtain 
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where  m  is  introduced  to  allow  for  the  application  of  Mark's  boundary 
condition  applicable  to  one- dimensional  situations  as  discussed  in 
Section  4.2.  We  shall  now  specialize  the  boundary  conditions  (6.2)  to 
the  specific  problems  of  harmonic  and  impulsive  motion  of  a  planar  wall, 
and  proceed  to  obtain  solutions  for  these  problems. 
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Harmonic  Motion  of  a  Planar  Wall.  The  problem  of  the  propagation  of 
disturbances  generated  by  harmonic  oscillation  of  a  planar  vail  was 
considered  by  Vincenti  and  Baldwin  (1962)  with  the  exponential  approx¬ 
imation.  We  shall  now  reconsider  the  same  problem  on  the  basis  of  the 
present  differential  approximation.  In  order  to  compare  our  results 
with  those  of  Vincenti  and  Baldwin,  it  is  convenient  at  this  point  to 
change  to  their  variables.  For  this  purpose,  we  let 


e  , 
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ba  a. 


16  T;  16  b, 

a  =  3Bo 


(6.5) 


where  cd  is  the  radian  frequency  of  oscillation  and  b  is  the  constant 
used  in  the  exponential  approximation.  Thus  b  *  1.562  in  Vincenti 
and  Baldwin's  work,  while  b  =  V5  for  the  P^-approximation  as  shown 
in  Section  4.3. 

The  boundary  conditions  (6.2)  for  harmonic  motion  can  be 

written 


u'(0,t)  =  (0, t)  =  Re[A  e^*1]  , 

°  *S 


(6.6a) 


(6.6b) 


where  A  and  B  are  dimensionless  complex  constants  assumed  to  be 
given. 

We  assume  the  solution  in  the  form 
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(6.7) 


FT 

0({,t)  -  -£■  Re[H(l)  etot]  , 

where  H(£)  is  a  function  to  be  determined.  Since  Eq.  (6.7)  is  periodic 
in  t,  substitution  of  this  equation  into  Eq.  (6.1)  leads  to  a  fourth- 
order,  constant-coefficient,  linear  ordinary  differential  equation  for 
H(£).  The  solution  of  this  equation  is  of  the  form 

4  c  A 

H(0  =  £  C.  e  J  .  (6.8) 

J-l  J 

Substituting  Eqs.  (6.7)  and  (6.8)  into  Eq.  (6.1),  we  obtain 

(Cj  -  62)  (1+Cj)  -  iKpr'1(/^Cj)  Cj  =  0  .  (6.9) 

Since  Eq.  (6.9)  is  quadratic  in  c^,  half  of  the  roots  will  have  a 
positive  real  part.  In  view  of  the  boundary  condition  (6.3),  these 
roots  are  inadmissible.  If  c^  and  c^  are  the  two  roots  with  pos¬ 
itive  real  part,  and  in  Eq.  (6.8)  must  be  taken  to  be  zero.  With 

these  considerations  and  imposing  boundary  condition  (6.6a),  we  find 

2 

2  c .  C .  =  A  .  (6.10) 

J=1  J  J 

Applying  boundary  condition  (6.4b)  and  utilizing  Eqs.  (6.6b)  and  (6.9), 
we  obtain 
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2  /Nr  ♦  =?)f  &  -  p  )  <c*  -  eV1  c.  -  b.  (6.iib) 

jti  J  Vyj  /  J  i 


If  Jtork's  boundary  condition  is  used,  Eq.  (6.11a)  with  m  *  V3  becomes 


S  r'1(r  +  c*)  3((Hc,)'1  C,  =  B  . 
J  J  J 


(6.11b) 


Equations  (6-9) ,  (6.10),  and  (6.11b)  are  sufficient  to  determine  all 
the  constants  in  the  solution  (6.7)*  As  one  would  anticipate  from  the 
analysis  in  Section  4.3,  these  equations  are  identical  to  Qqs.  (60), 
(6l),  and  (62)  in  Vincenti  and  Baldwin's  paper  (1962).  The  solution 
can  therefore  proceed  as  in  that  paper. 

Impulsive  motion  of  a  planar  wall.  The  solution  for  impulsive  motion 
of  a  planar  wall  has  been  obtained  previously  by  Baldwin  (1962)  and  by 
Lick  (1964),  again  using  the  exponential  approximation.  As  a  second 
example, we  shall  reconsider  here  this  problem  on  the  basis  of  the  differ¬ 
ential  approximation.  The  results  will  be  compared  with  those  of  Lick. 
It  should  be  noted  that,  although  Lick  solved  the  same  differential 
equation  (6.1),  his  boundary  condition  corresponding  to  Eq  (6.4b)  is 
in  integral  form. 

Following  Lick,  we  specialize  boundary  conditions  (6.2)  to 


t  <  0, 


u'(0,t) 
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(6.12a) 


Tw'(t)  =  0, 


(6.12b) 
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t  >  0 


i 


u'(0,t) 


(0,t)  =  0  , 


(6.12c) 


V(t) 


6 

R 


(6. 12d) 


The  initial  conditions  are 


0(x,O)  =  (x,0)  =  ^  (x,0)  =  0  .  (6.13) 

dt  at-2 

Let  the  Laplace  transform  of  0  be 

0(x,  p )  *  /“  e‘pt  0(x,t)  dt  .  (6.14) 

0 

Applying  the  Laplace  transformation  to  Eq.  (6.1)  and  imposing  initial 
conditions  (6.13),  we  obtain  a  fourth-order,  constant-coefficient, 
linear  ordinary  differential  equation  for  ?  The  solution  of  this 
equation  is  of  the  form 


#(x,  P)  =  2  Ai  e  ^  ’  (6.15) 

J=1  J 


Substituting  Eq.  (6.15)  into  the  ordinary  differential  equation,  we 
obtain 
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Equation  (6.16)  Is  quadratic  In  yy  Using  the  same  argument  as  In 
the  preceding  section,  we  retain  only  the  roots  with  negative  real  part. 
Let  these  two  roots  be  denoted  by  y ^  and  y Solving  the  algebraic 
equation  (6.16),  we  find 
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In  the  expressions  (6.17b),  b  =  3/2  was  used  by  Lick  in  his  exponential 
approximation  whereas  b  =  Y3  is  for  the  present  differential  approx¬ 
imation. 

Imposing  boundary  condition  (6.12)  on  (6.15),  we  obtain 


2 

S  k\  y  \  =  p/p  • 

j=i  J  J 


(6.18) 


Applying  the  Laplace  transform  to  boundary  conditions  (6.5b)  am  (6.13c), 
and  imposing  these  conditions  on  Eq.  (6.15),  we  also  obtain 


6  = 


i  2'  2, 
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(6.19a) 


If  Mark’s  boundary  condition  is  used,  Eq.  (6.19a)  with  m  =  V3  becomes 


2  A  (p2-a 2y2) 
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It  follows  from  Eqs.  (6.l8)  and.  (6.19b)  that 
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C1  S  bl(p2  *  ftT  7l)/(?l  *  V' 
C2  s  b^(pt  -  a2  72)/^2  +  bl^* 


(6.21) 


(6.22) 


Equations  (6.15)  and  (6.17)  with  expressions  (6.20)  through  (6.22)  are 
equivalent  to  Eq.  (3-2)  with  expressions  (3*4)  through  (3*7)  of  Lick's 
paper  (1964).  Thus  Lick's  integral  boundary  condition  is  identical  to 
the  differential  boundary  condition  (6.4)  with  m  =  V3  if  b  in  his 
exponential  approximation  is  suitably  chosen.  Again,  this  conclusion 
is  to  be  anticipated  from  the  analysis  in  Section  4.3* 
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6.2.  Stead  Flow  Over  a  Wavy  Wall. 

To  provide  a  simple  analytical  solution  for  multi- dimens iona  1 
radiating  flow,  we  now  consider  the  steady  two-dimensional  flow  over  a 
wavy  wall.  So  far  as  the  author  is  aware,  this  is  the  first  solution 
obtained  for  a  multi-dimensional  problem  with  no  restriction  on  the 
magnitude  of  the  absorption  coefficient  or  temperature  range. 

Ve  assume  that  the  flow  field  is  described  by  a  perturbation 
on  a  uniform  parallel  flow  with  velocity  Uw  in  the  x^-direction. 
Suppose  that  the  avy  wall  is  designated  by 


*2 


€  sin  2«  -y 


(6.23) 


where  c  is  the  amplitude  and  l  the  wave  length  of  the  wall.  To 
simplify  the  analysis,  the  wall  is  assumed  to  be  black  and  to  be  main¬ 
tained  at  a  temperature  that  is  a  small  departure  from  that  of  the  gas 
at  infinity.  The  governing  equation  for  this  case  is  given  by  Bq.  (5.16) 
specialized  to  two-dimensional  flow  as  follows: 
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(6.24) 


where 
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Within  the  framework  of  the  linearized  theory,  the  boundary 
conditions  at  the  wall  are 

,  dti  dXp  x, 

^(x^,  0)  *  (x^,0)  =  »  2nUe  cos  2n  -j-  ,  (6.25) 

dT'  2nT  x 

diT  ^Xl^  =  T~  D  cos  2n  T  »  (6.26) 

where  is  the  perturbation  wall  temperature  and  D  specifies  the 

amplitude  of  the  temperature  variation  on  the  wall.  The  fact  D  in 
general  be  complex  allows  for  a  possible  phase  shift  between  the  wall 
and  the  variation  of  wall  temperature.  The  factor  2 nT^/i  in  Eq.  (6.26) 
is  introduced  to  make  D  dimensionless.  At  infinity,  we  require  that 
the  disturbances  are  finite.  Thus  we  have 


0(x^,  »)  =  finite  quantity  .  (6.27) 


We  now  obtain  the  boundary  condition  (6.26)  in  terms  of  0. 
Specializing  Eq.  (5*l8b)  to  the  two-dimensional  case,  we  obtain  the 
perturbation  wall  temperature  in  terms  of  0  as 
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Imposing  the  boundary  condition  (6.26)  in  Eq.  (6.28),  we  have 
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The  periodic  nature  of  the  boundary  conditions  suggests 
that  the  solution  is  of  the  form 
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(6.30) 


where  the  factor  U  i/2n  is  introduced  to  make  H  dimensionless.  The 
dimensionless  independent  variables  x  and  y  are  related  to  the 
origninal  variables  by 

2nx,  2jcx^ 

x  =  — j —  and  y  =  — j — .  (6.31) 


Substitution  of  Eq.  (6.30)  into  Eq.  (6.24)  leads  to  a  fourth-order, 
constant-coefficient  ordinary  differential  equation  in  H(y).  The 
solution  is  of  the  form 


4  c  y 

H(y)  =  2  Ai  e  ^  >  (6.32) 

J=l  J 
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where  the  A^'e  and  c^'s  are  complex  quantities.  To  find  the  c^'s, 
we  substitute  Eqs.  (6. 30)  and  (6.32)  into  Eq.  (6.24)  to  obtain 


0Bu  /  2 
nBo  CJ 


+  «Vb8> 


H#+i]) 


(6.33) 


where  tc  -  1 

called  the  Bouguer  number.  Equation  (6. 33)  is  a  fourth-order  equation 
for  all  finite,  non-zero  values  of  Bo  and  Bu.  It  remains  fourth- 
order  when  Bo  =  0,  Bo  -*<»,  and  Bu  =  0.  When  Bu  approaches  infinity, 
however,  Eq.  (6.33)  reduces  to  a  second-order  equation.  The  roots  of 
c.  for  the  limiting  cases  of  the  parameters  Bo,  3u  and  M  are 

J  boo 

as  follows: 
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and  Bu  *  a  /.  The  quantity  Bu  is 
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and  Cj  =  1  , 

and  Cj  =  1  +  3(Bu/2«)2, 

and  c2  =  1  ,  (6.34a) 
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It  is  also  of  interest  to  note  that  for  y 
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In  the  general  case,  since  Eq.  (6.33)  is  a  quadratic  equation  in  c^, 
half  of  the  roots  will  have  positive  real  part.  In  view  of  the  boundary 
condition  (6.3)*  these  roots  are  inadmissible.  The  roots  with  negative 
real  part  can  be  found  by  solving  Eq.  (6.33)  formally  to  obtain 
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(6.3Ub) 

The  Aj's  in  Eq.  (6.32)  can  be  found  by  imposing  the  boundary 
conditions  as  follows:  Substituting  Eqs.  (6.30)  and  (6.32)  into  boundary 
condition  (6.25),  we  have 
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Substituting  Eqs.  (6.50)  and  (6.52)  into  boundary  condition  (6.29)  and 
making  use  of  Eq.  (6.55)>  we  also  obtain 


(6.56) 


To  simplify  the  analysis  we  now  assume  D  =  0  (i.e.,  the  wall  temperature 
is  constant  and  is  the  same  as  that  at  infinity).  Solving  the  simulta¬ 
neous  Eqs.  (6.55)  and  (6.56)  with  D  *  0,.  we  find 
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It  is  convenient  to  introduce  the  real  and  imaginary  parts 
of  c4  such  that  c^  *  -(5^  +  iX^),  where  64  and  are  positive 
quantities.  The  velocity  potential  Eq.  (6. 50)  as  a  function  of  x. 


and  x  is  then  given  by 
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where  a.  and  b  are  the  real  and  imaginary  parts  of  A./(2nc//). 

J  J  J 

Equation  (6.38)  represents,  in  general,  two  systems  of  waves  with  differ¬ 
ent  amplitude,  damping  factor,  and  direction  of  propagation.  Bich  of 
these  waves,  with  amplitude  Ay  decays  exponentially  with  x ^  with  a 

damping  factor  5.  and  along  straight  lines  with  slope  (measured  from 

J 

the  positive  xy coordinate)  of 


(6.39) 


The  variation  of  6y  and  as  functions  of  the  parameters 
Bu,  and  l  3  will  be  examined  in  detail  below.  We  shall  see  that  the 
Ayterm  differs  only  slightly  from  the  classical  solution,  whereas  the 
A  -term  has  no  counterpart  in  classical  theory.  Following  Vincenti 
and  Baldwin  (1962),  we  refer  to  the  Ay  and  A  -terms  as  the  modified 
classical  wave  and  the  radiation- induced  wave,  respectively. 

The  disturbed  velocities  are  given  by 
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where 


Within  the  framework  of  the  linearized  theory,  the  pressure  coefficient 
on  the  wall  is 
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Equations  (6.40a)  and  (6.4l)  show  that  both  the  horizontal  disturbed 
velocity  and  the  pressure  distribution  at  the  wall  are  out  of  phase  with 
the  wall,  while  Eq.  (6.40b)  shows  that  the  vertical  disturbed  velocity 
is  in  phase. 

The  drag  coefficient  per  unit  wave  length  of  the  wall  can  be 

found  from 
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(6.45c) 


(6.45d) 


and  for  supersonic  flow 


(6.45e) 


(6.45f) 


(6.45g) 


(6.45h) 


Equations  (6.44)  >ud  (6.45)  are  identical  to  Ackeret's  classical 
solution  (1928)  for  flow  over  a  wavy  wall.  For  the  limiting  case  of 
an  infinitely  hot  gas  (Bo  =  0),  the  results  are  the  same  as  given  by 
Bqs.  (6.44)  and  (6.45)  but  with  the  isentropic  Mach  number  replaced  by 
the  isothermal  Mach  number.  For  this  problem  we  thus  see  that  the 
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radiation- Induced  wave  disappears  In  the  four  limiting  cases.  In  these 
limiting  cases,  the  decay  of  the  disturbance  velocities  is  zero  for 
supersonic  flow.  The  drag  coefficient  and  the  clockwise  rotation  (from 
the  positive  x^- coordinate)  of  the  lines  along  which  the  exponential 
decay  of  the  disturbance  velocities  takes  place  are  both  zero  for  sub¬ 
sonic  flow. 

Figure  2  through  Fig.  9  show  the  variation  of  8^,  X^,  Aj 
and  C.  as  functions  of  the  parameters  Mc  ,  Bu,  and  Bo  for  7  -  7/5* 

OB 

Figure  2  shows  the  variation  of  8^  and  X^  over  the  range  of  Mach 
number  Mc  from  0  to  2.0  for  specific  values  of  Bu  and  Bo. 

The  solid  lines  are  values  of  8^  and  X^  for  the  three  limiting  cases 
of  a  completely  cold  gas,  a  transparent  gas,  and  an  opaque  gas.  These 
are  also  corresponding  to  Ackeret's  classical  solution.  For  inter¬ 
mediate  values  of  Bu  and  Bo,  the  damping  factor  8^  is  never  zero 
over  any  part  of  the  Mach-number  range.  As  the  Mach  number  increases, 

8.  at  first  decreases  slowly  from  the  value  of  1.0  at  M„  «  0 
1 

(incompressible  flow),  decreases  rapidly  in  the  near-sonic  region,  and 
then  varies  only  slightly  in  the  supersonic  region.  The  value  of  X^, 
which  is  a  measure  of  the  clockwise  rotation  of  the  lines  along  which 
the  exponential  decay  takes  place,  is  zero  for  incompressible  flow.  Its 
value  is  small  in  the  subsonic  region  for  all  values  of  the  Bouguer  and 
Boltzmann  numbers.  It  begins  to  grow  rapidly  in  the  near-sonic  region 


and  continues  to  grow  in  the  supersonic  region. 

Figures  2a  and  2b  show  the  variation  of  8^  and  X^  versus 

Mc  at  discrete  values  of  Bo  for  the  thin-  and  thick-gas  respectively 
00 

As  the  Boltzmann  number  increases,  the  curves  of  8^  and  X^  shift  to 
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the  right  toward  the  classical  curves.  By  comparing  these  two  figures, 
we  note  that  a  thin  gas  (Bu  =01)  with  Bo  =  1.0  behaves  essentially 
like  a  completely  cold  gas;  whereas  a  thick  gas  (Bu  =  10.0)  with  the 
same  Boltzmann  number  behaves  like  an  infinitely  hot  gas.  This  is  so 
because  for  a  gas  at  moderate  Boltzmann  numbers,  if  the  absorption 
coefficient  is  small  as  in  the  thin-gas  case,  the  effect  of  radiation 
is  relatively  unimportant.  On  the  other  hand,  if  its  absorption  coeffi¬ 
cient  is  large,  the  effect  of  radiation  becomes  important. 

Figures  2c  and  2d  show  the  variation  of  6^  and  versus 

Mc  at  discrete  values  of  Bu  for  a  hot  and  cold  gas  respectively. 

As  Bu  increases,  the  curves  for  5^  and  first  move  away  from 

the  classical  curve;  with  further  increase  in  Bu,  they  eventually  move 
back  toward  the  classical  curves.  In  the  limit  of  infinite  Bouguer 
number,  they  again  coincide  with  the  classical  curves.  By  comparing 
Figs.  2c  and  2d,  we  note  that  for  Bo  =  0.1,  except  for  the  cases  when 
the  Bouguer  number  is  extremely  snail  or  large,  the  gas  behaves  essentially 
like  an  infinitely  hot  gas  (cf.  Figs.  2a  and  2b);  whereas  for  Bo  =  10.0, 
the  gas  behaves  like  a  completely  cold  gas. 

Figures  5a  and  5b  show  the  variation  of  6^  and  X^  for  the 
complete  range  of  Bouguer  number  for  specific  values  of  the  Boltzmann 
number  and  the  Mach  number.  To  emphasize  the  boundary- layer- like  behavior 
of  the  results,  the  horizontal  scale  for  the  Bouguer  number  from  10.0 
to  »  has  been  arbitrarily  chosen.  Thus  the  results  in  this  range  of 
the  Bouguer  number  are  qualitative  in  the  horizontal  direction  while 
remaining  quantitative  vertically.  The  boundary- layer- like  behavior 
becomes  pronounced  only  for  certain  combinations  of  parameters,  namely. 
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for  a  hot  gas  with  extremely  smell  or  Large  Bouguer  number.  For  a  fixed 
pair  of  values  of  Mach  and  Boltzmann  numbers,  the  values  of  5^  and  A^ 
are  qualitatively  symmetric  with  respect  to  Bouguer  number. 

Figures  3c  and  3d  show  the  variation  of  5^  and  A^  for  the 
complete  range  of  Bo  for  specific  values  of  Bu  and  Me  .  The 
horizontal  scale  for  Bu  in  the  range  of  10.0  to  »  is  again  arbi¬ 
trarily  chosen.  For  a  fixed  pair  of  values  of  Me  and  Bu,  the  values 

of  5^  and  A^  are  asymmetric  with  respect  to  Boltzmann  number.  Again, 

the  boundary-layer-like  behavior  is  observed  for  small  Boltzmann  number 
with  either  very  small  or  very  large  Bouguer  numbers. 

Figure  4  shows  the  variation  of  B2  and  Ag  over  the  Mach- 
number  range  of  0  to  2.0  for  specific  values  of  Boltzmann  number 
and  Bouguer  number.  We  see  that  for  the  radiation- induced  wave,  both 
6  and  Ap  are  practically  independent  of  Mach  number.  For  all  finite, 

non- zero  values  of  Bo  and  Bu,  the  values  of  &2  and  Ag  never 

become  zero.  Thus  the  radiation- induced  wave  is  always  damped,  and  the 
clockwise  rotation  of  the  lines  along  which  the  exponential  decay  takes 
place  is  always  non- zero  for  all  Mach  numbers. 

Since  the  values  of  Bp  and  Ap  are  practically  independent 

of  Mach  number,  a  plot  of  these  quantities  versus  Bu  at  specific  values 

of  Bo  (or  versus  Bo  at  specific  values  of  Bu)  for  a  specific  Mach 

number  will  be  qualitatively  representative  of  all  Mach  numbers.  Figures 

5a  and  5b  are  such  plots  for  M_  =2.0.  Figure  5*  shows  that  as  Bu 

s» 

varies  from  0  to  ®,  the  damping  factor  Bp  increases  monotonically 
from  1  to  »  while  Ap  increases  monotonically  from  0  to  *. 

Figure  5b  shows  that  as  Bo  varies  from  0  to  »,  the  damping  increases 
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slowly  and  then  increases  or  decreases  to  a  certain  asymptotic  value 

depending  on  the  value  of  Bu.  At  the  same  time,  V  increases  from 

zero  to  a  maximum  and  then  decreases  to  zero.  The  variation  of  \ 

with  respect  to  Bo  is  thus  roughly  symmetric. 

The  amplitudes  of  the  modified  classical  wave  and  the  radiation- 

induced  wave  depend  on  the  boundary  conditions.  Figure  6  shows  the 

variation  of  the  amplitudes  of  the  two  waves  for  the  Mach-number  range 

from  0  to  2.0  for  specific  values  of  Bu  and  Bo.  It  is  seen  that 

the  modified  classical  wave  predominates  at  all  values  of  M0  ,  Bu  and 

00 

Bo.  This  is  only  true  for  this  particular  problem,  and  would  not  hold 

for  other  situations  (for  example,  for  a  flat  wall  with  a  periodic 

temperature  distribution).  For  the  limiting  cases  of  a  transparent  gas 

(Bu  -»  0),  an  opaque  gas  (Bu  -»«),  and  a  completely  cold  gas  (Bo  ->»),  the 

amplitude  of  the  radiation- induced  wave  becomes  zero,  and  the  amplitude 

of  the  modified  classical  wave  is  identical  to  the  classical  solution 

represented  by  the  solid  lines.  It  is  seen  that  in  the  classical  solution, 

the  amplitude  becomes  infinite  at  M  =1.  For  the  limiting  case  of 

00 

an  infinitely  hot  gas  (Bo  =  0),  the  amplitude  of  the  radiation- induced 

wave  is  again  zero,  and  the  amplitude  of  the  modified  classical  wave 

becomes  infinite  at  =1.  For  finite,  non- zero  values  of  Bo  and 

00 

Bu,  the  amplitude  of  the  modified  classical  wave  never  becomes  infinite. 

For  fixed  values  of  Bu  and  Bo,  with  increasing  Mach  number  the 
amplitude  of  the  modified  classical  wave  increases  slowly  until  the 
near-sonic  region  is  reached,  where  it  rises  sharply  to  a  maximum  value 
and  then  declines.  In  the  supersonic  region,  |a^ |/(2ne//)  decreases 
slowly  with  increasing  Mg  . 
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Figures  7a  and  7b  are  plots  of  |A.|/(2ne/i)  versus  Bu  for 

J 

specific  values  of  Bo  and  Mc  .  As  Bu  varies  from  0  to  «*, 

O 

00 

<A^(2jte//)  for  subsonic  flow  increases  from  a  certain  value  to  a 
maximum  and  then  declines  back  to  the  original  value.  For  supersonic 
flow,  |A^  |/(2 ne/Z)  decreases  from  a  certain  value  to  a  minimum  and 
then  increases  back  to  the  original  value.  The  boundary- layer- like 
behavior  is  again  observed  for  a  hot  gas  with  extremely  small  or  large 
Bouguer  number. 

Figures  7c  and  7d  are  plots  of  |A  |/(2jrc/i)  versus  Bo  at 

J 

discrete  values  of  Bu  for  subsonic  and  supersonic  flows  respectively. 
As  Bo  varies  from  0  to  the  amplitude  of  the  modified  classical 
wave  decreases  for  subsonic  flow  and  increases  for  supersonic  flow. 

The  rate  of  decrease  or  increase  is  rapid  for  a  gas  with  extremely 
small  or  large  Bouguer  number. 

Figure  8  6hows  that  at  all  non- zero,  fihite,  Bouguer,  Boltzmann 
and  Mach  numbers,  the  drag  coefficient  never  becomes  zero  or  infinite. 
For  subsonic  flow,  it  is  small  for  all  values  of  Bu  and  Bo.  In  the 
near-sonic  region,  it  rises  sharply  to  a  maximum,  and  then  declines 
sharply  into  the  supersonic  region  where  it  continues  to  decline  slowly. 

Figures  9a  and  9b  show  the  variation  of  the  drag  parameter 

2 

Ci/(2n€/ir  versus  Bu  at  discrete  values  of  Bo  for  subsonic  and 
supersonic  flows  respectively.  Again  the  boundary- layer- like  behavior 
is  observed  for  a  hot  gas  with  extremely  small  or  large  Bouguer  number. 

Figures  9c  and  9d  show  the  variation  of  C^/  (2ne/  £)‘  as  a 
function  of  Bo  at  discrete  values  of  Bu  for  subsonic  and  supersonic 
flows  respectively.  With  increasing  Bo,  this  parameter  increases  from 
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zero  to  a  maximum  and  decreases  back  to  zero  for  subsonic  flow.  For 
supersonic  flow,  it  increases  monotonically  from  one  non-zero  value  to 
another.  The  rate  of  changes  of  the  parameter  is  rapid  for  a  hot  gas 
with  extremely  small  or  large  Bouguer  number. 


VII.  CONCLUDING  REMARKS 

The  spherical-harmonic  method  has  been  used  by  nuclear  physicists 
extensively  in  connection  with  the  design  of  nuclear  reactors  in  the 
past  decade.  It  has  been  shown  that  the  first  approximation  of  the 
spherical-harmonic  method  is  highly  accurate  for  such  neutron- transport 
problems. 

The  formulation  of  radiating  gasdynamics  with  the  spherical- 
harmonic  (differential)  approximation  appears  promising.  Although  we 
have  no  way,  as  yet,  of  assertaining  the  accuracy  of  the  approximation 
for  multi- dimensional  flow,  the  corresponding  approximation  in  one¬ 
dimensional  problems  is  known  to  be  extremely  accurate.  This  is  shown 
by  the  works  of  Pearson  (1964)  and  Ferziger  (1965).  For  the  problem 
of  the  radiation-resisted  shock  wave,  the  exact  numerical  solution  of 
Pearson  compared  extremely  well  with  the  solution  by  Heaslet  and  Baldwin 
(1963)  on  the  basis  of  the  exponential  approximation,  which  in  turn  is 
equivalent  to  the  differential  approximation.  For  the  problem  of 
radiative  transfer  between  parallel  plates,  Ferziger  shows  that  the 
differential  approximation  is  extremely  accurate  as  compared  with  an 
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exact  solution.  Furthermore,  the  present  formulation  recovers  the 
correct  thin-  and  thick-gas  approximations  in  the  appropriate  asymptotic 
situations.  All  of  these  facts  indicate  that  this  approach  is  at  least 
qualitatively  correct. 

The  grey-gas  approximation  is  obviously  unrealistic  from  a 
physical  point  of  view.  For  the  linearized  theory,  where  temperatures 
and  other  physical  properties  do  not  change  appreciably,  this  assumption 
may  not,  however,  be  critical. 

From  the  mathematical  point  of  view,  the  assumption  of  non¬ 
scattering  is  not  essential  to  the  differential  approximation.  Had  we 
taken  into  consideration  the  effect  of  scattering,  the  coefficient  in 
the  linearized  equations  would  be  modified  accordingly. 

Since  radiative  transfer  is  an  irreversible  process,  we  expect 
that  the  effects  of  radiative  nonequilibrium  might  be  similar  in  some 
ways  to  those  of  chemical  or  vibrational  nonequilibrium.  The  solution 
for  radiative  flow  over  a  wavy  wall  does,  in  fact,  show  the  occurrence 
of  pressure  drag  at  subsonic  speeds  and  a  smoothing  of  the  transition 
from  subsonic  to  supersonic  speeds.  In  these  respects,  the  results  are 
indeed  similar  to  the  flow  over  a  wavy  wall  with  chemical  or  vibrational 
nonequilibrium  (Vincenti  (l*?*^))* 
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APPENDIX  A 


RECURRENCE  AND  ORTHOGONALITY  RELATIONS 
FOR  THE  SPHERICAL  HARMONICS 


In  this  section,  we  shall  cite  without  proof  the  well-known 
formulas  and  relations  of  the  spherical  harmonics.  For  a  detailed 
discussion  of  spherical  harmonics,  the  reader  is  referred  to  standard 
texts  on  advanced  mathematics  (see  for  example,  Sommerfeld  (1964)  or 
Irving  and  Mullineaux  (1959))* 

The  spherical  harmonic  Y^(fl)  is  related  to  the  associated 
Legendre  function  P^cos  0)  by 

Y^(Q)  *  C”  ^  p”(w)  ,  (A-l) 


where  the  upper  Index  ffl  la  positive,  p  =  cob  0,  and  -  ^TTTmTT^^' 
If  the  upper  index  is  negative,  it  is  related  to  the  complex  conjugate 

Y^(fl)  by 


Y^m(fi)  =  (-1)”  c"  e'^  P™(u)  =  (-l)m  Y^(n).  (A-2) 


It  follows  from  (A-l)  and  (A-2)  that 


Y^(n)  .  T^(n) 


y^(n)  +  (-i)“  ./(a) 


(sin  0)m 


cos  P  <p)  , 

du 

(A- 3a) 
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and 


y^(5)  -  5*2(5)  y*2(5)  -  (-1)“  y;m(5)  .  d® 

/  /  1  i -  .  (sin  «)“  sin  up  5-  P,(n)  . 

m  i 


2C 


du 


(A-3b) 


In  view  of  Eqs.(A-l),  (A-2),  and  (A-3)»  the  direction  cosines  given  by 
Eq.  (2.18)  are,  in  terms  of  the  spherical  harmonics. 


i2  *  cos  e  *  yj(n)  , 


/.  *  sin  6  sin  0 


yj(fi)  -  ^(a) 

~ yi 


Y^(n)  +  Y"L(n) 

vi 


=  sin  0  cos  0  * 


Y^(Q )  4  ?J(5)  y*(5)  -  y"1^) 


V2 


V2 


The  addition  theorem  for  spherical  harmonics  is 


P/(ft  •  0')  =  S  Y^)  • 


nui- 1 


(A-U) 

(A- 4b) 


(A-4c) 


(A-5) 


The  orthogonal  property  and  recurrence  relations  for  the  spherical 
harmonics  are  as  follows: 


/  yjt'(n)  Y^(n)  an 
a 


(-i)m  /  y’“'(o)  y?(5)  dfi 

a  1  1 


4n 

T27TTT 


t 


(A-6) 
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APPENDIX  B 

DETAILED  DERIVATION  OF  EQUATION  (U.4) 

The  radiation-transport  equation  for  a  grey  gas  in  quasi¬ 
equilibrium  is 

^  ||  +  fl  *  grad  I  =  -(a+Tj)  I  +  2£!L  +  n  /  Kdfj*  .  (B-l) 

Following  the  procedure  used  in  neutron- transport  theory,  we  first 
expand  the  radiation  intensity  in  a  series  of  spherical  harmonics  as 
follows: 
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I(r,  ii,  t)  = 


2  2  t)  /I(o)  > 

1=0  HU-/ 


(B-2) 


where  the  A°(r,t)'s  are  functions  to  be  found  and  the  Y^(ft)'8  are 
spherical  harmonics  For  reasons  explained  in  Section  4.1,  the  scattering 
function  can  be  expanded  in  a  series  in  terms  of  the  Legendre  polynomials 
P.(cos  9  ).  Utilizing  the  addition  theorem  (A-5)>  we  can  therefore  write 

L  O 


*(•*„.  fi,)  =  2  Ki  PA0>  *  2  *t  Tt(S  •  s">  -  2  E 

/=0  t=0  1=0  m=-i 

(B-3) 


Substituting  Eqo.  (B-2)  and  (B-3)  into  Eq.  (B-l),  we  have 


<Ht  +fi  •  v  +  “  +  2  E  y>)  -  ^ 

1=0  m a-/ 


00  / 


»  / 


T)  /  2  E  A®(r,t)  An*)  2  E  *,  Y>)  Y^(ft')  dfi' 

fl  »  1=0  m =-/  Z  *  1=0  m  =-l  11  1 


(B-4) 


Using  the  orthogonality  relation  (A-6),  we  can  simplify  the  integral 
term  in  Eq.  (B-4)  to 


00  l 


n  2  E  A®(f,t)  2  */  E  ^^)  /  ??(«’)  «' 

i=0  m=-/  /— 0  m--/  fl ' 


i  .2.  E.  *>)  (jfe)  • 


1=0  m=-i 


(B-5) 
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With  this  result,  Eq.  (B-4)  can  be  rearranged  to  yield 


V  +  (a+Tj) 


‘/‘•"I 
1  2/+1J 


(B-6) 


To  eliminate  the  spherical  harmonics  in  Eq.  (B-6),  we  multiply 
Eq.  (B-6)  by  Y^(ii)  and  then  integrate  over  the  whole  range  of  solid 
angle.  The  first,  third,  fourth,  and  fifth  terms  in  Eq.  (B-6)  immediately 
yield 


iK 

c  5t 


(a+rj) 


Am  - 


X2I+1) 


Am  - 
A/ 


aaT 

n 


01  Om 


(B-7) 


The  second  term  in  Eq.  (B-6)  needs  more  attention.  We  first 

note  that 


a 


sin  6  sin  0 


t  cub  a 


dx~ 


dx^ 


(B-8) 


If  the  second  term  in  Eq.  (B-6)  is  multiplied  by  Y^(ft)  and  then 
integrated  over  the  whole  range  of  solid  angle,  it  becomes,  after  the 
substitution  of  Eq.  (B-8), 


l 

E 

m =-/ 


/  Y'J(fl)  Y^!(fl) 


sin  6  sin 


(B-9) 


+  cos  6 


+  sin  9  cos 
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Consider  now  the  first  term  in  Eq.  (B-9)*  With  the  aid  of  Eqs.  (A-l) 
and  (A-2-),  and  after  replacing  sin  0  by  its  exponential  form,  we  have 


dA® 

/ 

^  £ 


x[/ 


sin  9  e 


i(m-m'+'')0. 


_ir ,  . 

P,,(p) 


dft  +  /  sin  6 


i(m-m'-l)0  m,  .Tm',  ... 
e  '*yFi(p)P/,(M)df< 


(B-10) 


Making  use  of  Eq.  (A-9)  and  then  utilizing  Eq.  (A-7),  we  obtain  the  first 
term  in  Eq.  (B-10)  as 


-  i2rr 


[(/+m)(/+m-l)]1/t 

(2i-l)(2/+l) 


dA 


m-1 

£+1 


Similarly,  using  Eq.  (A-10)  and  Eq.  (A-7),  we  obtain  the  second  term  in 
Eq.  (B-10)  as 


-  i2n 


[(i-m)U-m-l)]1 

72£-l7  (2/+1) 


The  remaining  terms  in  Eq.  (B-9)  can  be  handled  in  a  similar  manner. 
Collecting  these  results  and  those  of  Eq.  (B-7)>  we  obtain  Eq.  (4.4). 


Figure  3.  Values  of  6,  and  X_  as  functions  of  Bu  and  Bo  at  subsonic  and  supersonic  speeds 


0.8J 


and 


Tirurm  6.  %1un  of  |A, |/(2«c/«) 


Figure  7.  Values  of  *A  ! /(2ne//)  as  functions  of  Bu  and  Bo  at  subsonic  and  supersonic  speeds 
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